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Objective 

The  long-term  objective  of  this  work  is  to  develop  a  quantitative  understanding  of  the  fluctuations  in  phase 
and  amplitude  of  seismic  wave  propagation  due  to  kilometer-scale  variations  in  wave  velocity  within  the 
earth.  The  theoretical  approach  is  to  describe  these  variations  in  a  statistical  way:  in  particular  to  consider 
the  variations  as  represented  by  a  spectrum  that  depends  on  depth,  and  may  depend  on  geographical  loca¬ 
tion.  Data  that  are  relevant  to  this  approach  include  wave-forms  with  frequency  content  above  one  Hertz 
received  on  seismic  arrays  or  on  world-wide  networks.  Relatively  high-frequency  data  is  desirable  because 
the  ability  to  discriminate  small  structure  is  dependent  on  the  wave  having  relatively  short  wavelength. 
Seismic  arrays  whose  elements  are  spaced  in  the  kilometers  to  tens-of-kilometers  regime  provide  analysing 
power  in  that  regime.  Arrays  with  larger  spacing,  such  as  world-wide  networks,  can  still  probe  small 
scales  if  the  various  available  sources  have  separations  in  the  above  range:  this  can  occur  for  earthquakes 
in  active  regions,  or  for  nuclear  explosions  distributed  within  test  sites. 
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A  realistic  understanding  of  the  small-scale  structure  in  the  earth  is  important  to  fundamental  geophysics 
because  it  affects  our  understanding  of  the  fundamental  dynamical  processes  in  the  earth.  Mantle  convec¬ 
tion,  chemical  differentiation,  fluid  permeation,  subduction-zone  dynamics,  and  crack  formation  all  will 
have  their  effects  on  small-scale  structure,  so  that  creation  of  more  sophisticated  theories  of  these  processes 
will  influence,  and  will  be  influenced  by,  our  understanding  of  small-scale  structure. 

It  is  not  likely  that  the  seismological  community  will  ever  have  a  complete  map  of  inhomogeneities  in  the 
earth  down  to  kilometer  scales.  Therefore  we  will  not  be  able  to  completely  predict  travel-time  and  ampli¬ 
tude  fluctuations  for  a  source-receiver  geometry  that  is  even  a  few  kilometers  different  from  previously 
measured  situations.  If  we  have  a  realistic  statistical  picture  of  the  small-scale  structure  in  wave-speed 
within  the  earth,  then  the  theory  of  wave  propagation  through  random  media  (WPRM)  can  be  used  to 
predict  the  scale  and  strength  of  travel-time  and  amplitude  fluctuations  due  to  earth  structure.  This  infor¬ 
mation  can  then  be  used  to  calculate  the  accuracy  of  yield  estimates  and  detection  thresholds  based  on 
seismic  information  from  an  arbitrary  array  of  seismometers,  with  a  priori  knowledge  of  results  from 
nearby  explosions  or  earthquakes.  Furthermore,  this  knowledge  can  be  used  to  design  arrays  in  an  optimal 
fashion  to  deal  with  random  earth  structure. 

Conferences  attended 


The  principal  investigator  attended  the  following  meetings  during  the  contract  period: 


Acoustical  Society  of  America  Meeting  in  Seattle,  May,  1988.  A  talk  on  supercomputer  simulation  of 
wave  propagation  was  given. 


Massachusetts  Institute  of  Technology,  talk,  June,  1988.  A  talk  on  supercomputer  numerical  simulation, 
and  applications  to  seismology  was  given. 

Jason  summer  study,  July  1-31,  1988.  The  principal  investigator  led  a  sub-group  investigation  of  regional 
wave  propagation  analysis. 


Workshop  on  Wave  Propagation  in  Random  Media,  in  Tallin,  Estonia,  USSR,  September  19-23,  1988.  A 
talk  on  supercomputer  simulation  of  wave  propagation  was  given. 


n  For 
\*I 


Institute  of  Atmospheric  Sciences,  Moscow,  September  26, 1988.  Discussions  with  Academician  Obukhov 

_  •  Mon., 

and  Professor  V.I.  Tatarskii  were  carried  out.  Their  interests  include  general  wave  propagation  in  random  _ 


media. 


Lebedev  Institute  of  General  Physics,  Moscow,  September  26,  1988. 
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Discussions  with  Drs.  Bunkin, 
mat  : 


1  on/ 

titr  Socles 
I  and/or 
Cpectal 


□  □ 


Shishov,  and  Kravtsov  were  carried  out.  Their  interests  include  the  use  of  the  theory  of  wave  propagation 
in  random  media  in  ocean-acoustic  applications. 

Dr.  Jan  Martin  attended  a  meeting  of  1UGG  (the  International  Union  of  Geodesy  and  Geophysics)  in  Van¬ 
couver,  Canada,  August  9-22, 1988,  where  he  gave  a  talk  and  discussed  the  use  of  the  theory  of  wave  pro¬ 
pagation  in  random  media  in  plasma  physics  applications. 
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Inhomogeneities  Near  the  Core-Mantle  Boundary  Inferred  from 
Short-Period  Scattered  PKP  Waves  Recorded 
at  the  Global  Digital  Seismograph  Network 

Klaus  Bataii.uk  and  Stam.ky  M.  Flattl 


Physics  Department  and  Institute  of  Tectonics.  University  of  California.  Santa  Cruz 


We  analyze  short-period  PKP  precursor  wave  trains  recorded  al  Global  Digital  Seismograph  Network  sia 
lions  in  the  distance  range  1 20°  <  A  <  140°  to  infer  the  nature  of  inhomogcneitics  near  the  core-mantle 
ooundary.  Travel  limes  and  particle  motions  are  consistent  with  predictions  of  single-scattered  waves  near 
the  CMB.  The  dominant  frequencies  are  around  1  llz.  The  regions  best  sampled  arc  beneath  Indonesia, 
North  America.  Central  Afnca,  and  South  America  Based  on  first-order  scattering  theory,  the  characteris¬ 
tics  of  the  wave  number  spectrum  of  the  structural  inhomogcneitics  arc  obtained  for  two  hypotheses: 
volume  inhomogeneiucs  and  topographic  irregularities,  for  the  range  of  wavelengths  sampled  by  these 
data  (between  10  and  70  km),  the  spectrum  of  inhomogcneitics  is  best  represented  by  a  power  law  of  index 
5.3  for  volume  inhomogencities  and  6.8  for  topographic  irregularities.  If  the  scattering  is  due  only  to 
volumclnc  inhomogcneitics,  we  estimate  the  strength  of  the  structural  variations  to  be  about  0.5%  in  P 
wave  velocity  for  a  200-km-thick  layer.  If  the  scattering  is  due  only  to  topographic  irregularities,  the 
topography  has  an  rms  height  of  280  m.  At  present  we  can  not  distinguish  between  these  two  types  of 
inhomogcneitics. 


1.  Introduction 

There  is  increasing  evidence  for  the  existence  of  hetero¬ 
geneities  near  the  core-mantle  boundary  (CMB)  which  span 
wavelengths  from  approximately  10  to  3000  km.  The  upper  limit 
of  these  wavelengths  is  intrinsically  constrained  by  the  size  of 
the  core,  and  the  lower  limit  is  observational!)'  constrained  by 
the  probing  wavelength  of  seismic  waves.  Although  the  evi¬ 
dence  for  heterogeneities  is  well  established,  there  is  still  consid¬ 
erable  debate  about  its  physical  nature.  Among  the  proposed 
models  for  the  nature  of  inhomogeneiucs  are  (1)  thermal  boun¬ 
dary  layer,  (2)  chemical  boundary  layer,  and  (3)  topographic 
irregularities.  Currently,  it  is  not  possible  to  rule  out  any  of  these 
models  from  seismological  observations  alone. 

There  is  agreement  that  the  temperature  difference  between 
the  lop  of  D"  and  the  CMB  is  al  least  1000°K  [Stacey  and 
Lopcr.  1983;  Jeanlnz  and  Richter,  1979|;  therefore  a  thermal 
boundary  layer  must  be  developed.  Yuen  and  Peltier  [198(1) 
suggested  that  the  thermal  boundary  layer  is  unstable  and  hence 
dynamically  active,  and  not  just  the  edge  of  a  convective  cell. 
Stacey  and  Loper  (1983|  suggested  a  mechanism  for  how  mantle 
plumes  originate  at  the  bottom  of  the  mantle,  providing  a  solu¬ 
tion  of  the  temperature  distribution  within  D"  near  the  plume. 
This  solution  is  consistent  with  observations  of  long-period 
core-diffracted  waves  [Doornbos  et  at.,  1986).  Attempts  at 
numerical  simulation  of  this  high-Rayleigh  number  convection 
situation  are  being  made  now  (e.g.,  Boss  and  Sacks  [  1 985  J  and 
Olson  et  al  1 1987a,  b ]  among  others),  as  well  as  laboratory 
experiments  [e.g..  Loper  and  McCartney,  1986|. 

Two  models  of  chemical  boundary  layers  arc  proposed.  First, 
due  to  the  large  density  contrast  between  the  core  and  the  man¬ 
tle,  materials  with  densities  smaller  than  the  core  and  larger  than 
the  mantle  would  remain  at  the  CMB  [Ringwood.  1 979 [ . 
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Second,  based  on  the  phase  diagram  of  FcO  at  high  pressures 
and  temperatures,  Kniltle  and  Jeanlnz  [  1986]  have  suggested 
that  the  metallic  iron  core  reads  chemically  with  the  oxygen  of 
mantle  silicates,  with  the  mantle  residuum  of  such  a  reaction 
corresponding  to  the  D“  layer. 

Among  the  seismological  evidence  for  the  presence  of  inho¬ 
mogeneities,  Ritzwoller  el  al.  [1986]  find  indications  of  large- 
scale  inhomogeneitics  in  the  core  or  lower  mantle  from 
anomalous  splitting  of  cigenfrcquencies  associated  with  PKP , 
SKS ,  and  PKIKP .  Poupinet  el  al.  [1983],  Morelli  and  Dz'tewon- 
ski  [1987],  and  Creager  and  Jordan  [1986]  also  provide  evi¬ 
dence  of  large-scale  inhomogeneitics  near  the  CMB  from  inver¬ 
sion  of  travel  times  of  PKP  (BC  and  DF)  phases  for  thousands 
of  earthquakes.  Different  models  for  the  structure  near  the  CMB 
have  been  obtained  from  studies  of  waveforms  and  apparent 
slownesses  of  long-period  core  diffracted  P  and  S  waves  |A/<m 
ander  and  Phinney.  1966;  Doornbos  and  Mondt,  1979;  Mula 
and  Muller,  1980]  and  S ,  SKS,  and  ScS  waves  [Mitchell  and 
Helmberger,  1973;  Lay  and  Helmberger,  1983;  Schlillenhardl. 
el  al.  1985),  some  with  smooth  variations,  others  with  sharp 
discontinuities.  Differences  in  the  proposed  models  for  the  struc¬ 
ture  of  D"  might  reflect  its  strong  lateral  variation.  The  effect  of 
lateral  variations  on  some  of  these  phases  is  still  debated  [Cor¬ 
mier,  1986J. 

Variability  of  small-scale  heterogeneities  can  be  tested  by 
studying  short-period  PKP  waves  scattered  at  different  locatioas 
worldwide.  Itaddon  and  Cleary  [1974],  King  et  al.  |1974], 
Husebye  el  al.  [1976],  and  others,  presented  clear  evidence  that 
the  short-period  wave  train  arriving  before  PKIKP  originates 
from  scattering  by  irregularities  near  the  CMB.  They  based  their 
evidence  on  studies  of  data  from  large  seismic  arrays.  Although 
these  authors  find  differences  between  regions  sampled  by  their 
data,  they  agree  that  a  few  percent  in  density  and  seismic  veloci¬ 
ties  in  the  lower  mantle  (D "),  or  a  rough  CMB  with  radial  varia¬ 
tions  up  to  a  few  hundred  meters  could  produce  the  observed 
energy  level  of  precursors  to  PKIKP  j Van  der  Berg  et  al. 
198()|.  Here  we  report  on  an  analysis  using  available  source- 
receiver  geometries  during  a  period  of  6.5  years  recorded  al  the 
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Fig.  1.  Paths  of  PKIKP  (thin  solid  curve),  direct  PKP  (thick 
solid  curve),  and  scattered  PKP  (thick  dashed  curve)  at  entrance 
and  exit  from  the  core  arriving  in  the  distance  range  120° 
<  A  <  150°. 


Global  Digital  SemogTaph  Network  (GDSN).  Our  aim  is  to 
determine  the  nature  of  inhomogeneities  for  each  region  sam¬ 
pled.  In  section  2  we  present  the  data  and  their  characteristics; 
in  section  3  we  compare  the  data  with  the  predictions  of  two 
scattering  models:  volumetric  inhomogene i ties  and  topographic 
irregularities;  and  in  section  4  we  discuss  their  relation  and  some 
geophysical  implications. 

2.  Data 

Precursors  to  PKIKP  in  the  distance  range  i20°  <  A  <  140° 
recorded  at  the  GDSN  between  1980  and  1986  were  analyzed. 
The  paths  of  PKIKP  and  scattered  PKP  arriving  as  precursors 
are  shown  in  Fig.  1.  A  rough  estimate  is  that  200  events/year 
with  magnitudes  mk  >  5.7  occur  in  the  world.  In  principle,  each 
of  these  events  should  have  at  least  one  GDSN  station  within  the 
appropriate  distance  range  for  recording  precursors.  Thus,  in 


principle,  a  total  of  1300  events  are  expected  during  the  6.5 
years  of  data  analyzed.  Judged  from  their  clear  precursor  signal, 
130  records  at  20  stations  from  100  earthquakes  were  selected 
among  available  source -receiver  paths. 

Using  data  from  large  seismic  arrays  it  is  possible  to  deter¬ 
mine  the  slowness  of  the  arriving  signal  and  distinguish  in  this 
way  precursors  (as  scattered  waves  near  the  CMB)  from  noise. 
With  a  single-station  method  this  is  not  possible.  We  distin 
guish  precursor  signals  from  noise  based  on  three  main  criteria 
First,  to  clearly  identify  the  PKIKP  arrival,  events  with  sharp 
PKIKP  onsets  were  selected.  In  general  this  is  the  case  and  only 
5%  of  the  total  possible  data  set  were  discarded  for  this  reason. 
Second,  we  compared  the  amplitude  and  frequency  content  of 
precursors  with  the  noise  background  from  more  than  25  s 
before  PKIKP .  Amplitude  of  precursors  above  twice  the  noise 
background  were  selected;  this  accounted  for  about  85%  of  the 
selected  set.  For  smaller  amplitude  ratios,  we  compared  the  fre¬ 
quency  content;  typically,  precursors  have  a  peak  between  1  and 
2  Hz  while  the  noise  has  predominantly  lower-frequency  con¬ 
tent;  this  criteria  accounted  for  about  5%  of  the  selected  set. 
Third,  the  particle  motion  was  examined  when  three-component 
records  were  available;  only  signals  with  sleep  angle  of 
incidence  were  selected,  which  accounted  for  another  10%  of 
the  selected  set.  This  third  criterion  is  particularly  useful  for 
small  distances  where  the  amplitude  of  precursors  could  be  as 
small  as  one  fiftieth  the  amplitude  of  PKIKP .  It  will  be  shown 
later  that  the  large  number  of  events  in  the  correct  angular 
region  that  were  not  selected  (because  they  showed  no  visible 
precursor)  are  also  consistent  with  our  model  of  CMB  scattering, 
because  our  model  predicts  their  precursor  amplitude  to  be 
below  the  noise  level.  A  more  detailed  description  of  the 
analysis  and  all  the  waveform  data  can  be  found  in  the  work  by 

Bataille  [1987). 

The  regions  of  the  CMB  that  are  sampled  by  the  data  are 
shown  in  Fig.  2.  The  path  between  Indonesia  and  North  America 
is  better  sampled  due  to  the  seismically  active  region  in 
Indonesia  and  the  good  station  coverage  in  North  America. 


Fig.  2.  Regions  of  the  CMB  sampled  by  the  data.  Dark  (light)  shaded  regions  represent  scattering  at  the  entrance 
(exit)  of  the  core.  Large  light  diamond  represent  receivers,  and  small  dark  diamonds  represent  sources. 
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Fig.  3.  An  earthquake  in  Indonesia  (9.4°  N,  7.7°  W)  recorded  at 
a  distance  of  53.6°  (SNZO)  showing  a  simple  P  wave  arrival, 
while  at  138.5°  (RSNY)  and  140°  (SCP)  a  clear  precursor  wave 
train  is  observed.  The  total  display  time  is  90  s. 

One  typical  example  of  the  precursor  wave  train  for  one  earth¬ 
quake  observed  at  different  stations  is  shown  in  Fig.  3.  The 
duration  of  the  precursor  and  its  uniform  amplitude  along  the 
wave  train  are  characteristics  that  can  be  explained  by  a  scatter¬ 
ing  mechanism.  For  the  event  shown  in  Fig.  3,  the  station 
RSNY  has  three  component  records  available  (Fig.  4(top)),  and 
one  can  see  from  its  particle  motion  diagram  (Fig.  4(bollom }} 
that  the  incidence  angles  of  precurr  or  and  of  PKIKP  do  not 
differ  by  more  than  10°.  This  small  difference  implies  that  pre¬ 
cursors  have  passed  through  the  CMB.  Here  we  do  not  question 
further  the  origin  of  precursors  and  assume  that  all  our  data 
selected  in  this  way  are  due  to  scattering  near  the  CMB. 

Variations  of  the  properties  of  inhomogeneities  from  region  to 
region  would  be  manifested  in  variations  of  properties  of  the 
precursor  wave  train  as  recorded  in  different  stations.  Two  pro¬ 
perties  of  the  precursor  wave  train  can  be  compared:  the  com¬ 
plexity  of  the  waveform,  and  its  amplitude  compared  to  PKIKP  . 
All  precursor  wave  trains  have  in  general  a  continuous  ampli¬ 
tude  as  a  function  of  arrival  time  (see,  for  example,  Figure  3); 
there  are  only  5%  of  these  events  with  peaks  localized  in  time. 
These  peaks  in  amplitude  might  correspond  to  strong  construc¬ 
tive  interference  among  all  scattered  waves  contributing  at  that 
particular  time,  or  to  a  stronger  inhomogeneity  at  a  particular 
patch  wiLhin  the  scattering  region.  With  a  single  station  it  is  not 
possible  to  resolve  these  two  alternatives.  We  filter  the  PKIKP 
and  precursor  waveforms  between  0.6  and  4  Hz  using  a  Butter- 
worth  filter  of  first  order.  We  do  not  find  a  systematic  difference 
in  the  frequency  content  of  the  precursor  signal  as  a  function  of 
epicentral  distance,  or  as  a  function  of  geographical  location. 
After  filtering  the  PKIKP  and  precursor  signal  at  the  predom¬ 
inant  frequency  of  PKIKP ,  we  compute  the  amplitude  of  precur¬ 
sor/  PKIKP ,  which  as  a  function  of  distance  is  shown  (for  the 
stations  with  the  most  data  points)  in  Fig.  5.  All  stations  have  a 
large  scatter  of  data,  in  particular  stations  ZOBO  and  BCAO. 
This  scatter  could  correspond  to  differences  in  the  source 
geometry,  source  spectrum,  or  path.  However,  using  the 
moment-tensor  solution  from  the  Harvard  catalog,  we  compared 
the  radiation  pattern  of  PKIKP  and  their  precursors,  and  varia- 


Fig.  4.  (Top)  Vertical  (Z),  radial  (R),  and  transverse  (T)  com¬ 
ponent  records  a'  RSNY  (A  =  138.5°).  Note  that  the  scales  are 
not  the  same.  (Bottom)  Particle  motion  for  four  time  windows: 
(a)  noise  preceding  the  precursor  wavetrain,  (b)  precursor 
wavetrain,  (c)  PKIKP  ,  and  (</)  PKIKP ’s  "coda." 


tions  in  general  were  negligible.  For  the  events  during  1983  that 
were  in  the  correct  distance  range,  but  that  we  did  not  select 
because  no  visible  precursors  were  observed,  we  plotted  the 
noise  level  during  1983  with  a  diamond  symbol  in  Fig.  5.  It  can 
be  seen  for  those  cases  that  the  noise  level  falls  above  the  aver¬ 
age  level  of  scattered  waves  and  therefore  is  consistent  with  our 
model  of  CMB  scattering.  There  are  two  exceptions  at  station 
ZOBO,  but  that  station  has  the  largest  scatter,  and  we  see  no 
inconsistency  with  our  model.  Although  the  data  are  limited,  we 
find  no  systematic  variation  of  the  amplitude  ratio  at  a  given  sta¬ 
tion  as  a  function  of  azimuth.  Also,  the  amplitude  ratios  at  a 
given  distance  are  equal  for  all  stations  within  the  scatter  of  the 
data.  We  have  found  no  model  of  regional  variation  in  the 
strength  of  inhomogeneities  over  the  surface  of  the  CMB  that 
reduces  the  scatter  in  the  data.  We  conclude  that  the  spectral 
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Fig.  5.  Amplitude  ratio  of  precursor /PKIKP  as  a  function  of  epicentral  distance.  The  asterisks  are  calculated  with 
events  in  which  a  precursor  is  observed.  For  events  in  which  no  precursor  is  identified,  the  diamond  upper  limits  are 
calculated  from  the  noise  level  during  the  precursor  arrival  time. 


strength  of  the  inhomogeneities  is  the  same  from  region  to 
region  to  within  a  factor  of  three.  Van  den  Berg  el  al.  [1980] 
measured  the  amplitudes  of  precursors  recorded  at  4  UKAEA 
arrays  and  NORSAR.  Their  observations  of  the  ratio  of  precur¬ 
sor  amplitude  to  PKIKP  at  NORSAR  and  three  of  the  UKAEA 
arrays  are  consistent  with  our  results  in  Fig.  5.  One  of  the  arrays, 
YKA,  shows  consistently  lower  ratios,  indicating  less  scattering. 
The  data  at  YKA  map  the  regions  beneath  the  S.  Sandwich 
Islands  (at  entrance)  and  Central  North  America  (at  exit),  and  all 
have  an  epicentral  distance  of  about  136°.  One  tentative 
interpretation  is  to  pinpoint  the  region  beneath  Central  North 
America  as  less  inhomogeneous  compared  to  the  other  regions 
sampled,  but  more  data  is  necessary  to  confirm  this  observation. 
In  the  next  section  we  specify  in  more  detail  the  implications 
that  our  observations  have  on  the  properties  of  inhomogeneitics. 

3.  Result  Analysis 

Our  interpretation  of  precursors  to  PKIKP  is  based  on  first- 
order  scattering  theory.  In  first-order  theory  each  point  of  a 


heterogeneous  medium  becomes  a  source  of  scattered  waves  of 
both  P  and  5  types.  Only  single-scattered  waves  are  taken  into 
account.  A  brief  derivation  of  the  formulas  used  below  is 
presented  in  the  appendix.  The  mean  square  amplitude  of  scat¬ 
tered  waves  can  be  expressed  as 

<lul2>  =  — t*f  A\x’)R\x')Q.(x')<t>(K)d1x'  (1) 

4  JCMS 

where  k  is  the  wavenumber  of  the  wave  field,  A  is  the  product 
of  the  geometric  spreading  coefficient  from  the  source  to  the 
scattering  point  and  the  geometric  spreading  coefficient  from  the 
scattering  point  to  the  receiver,  R  contains  all  the 
refleclion/transmission  coefficients  and  attenuation  factors, 
K  =  k(ri  -  ri')  where  ri  and  ri'  are  the  direction  of  the  incident 
and  scattered  wave  at  the  scattering  point,  <t>(K)  is  the  spectrum 
of  inhomogeneities,  ft  =  1  for  volume  inhomogeneities,  and  for 
topographic  irregularities  ft  has  a  more  complicated  expression 
(see  the  appendix  equation  (Al 8)).  Assuming  9  to  be  the 
scattering  angle,  we  have  for  volume  inhomogeneities 
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5600=  Mil  co  S0-^1L 
p  X  +  2  4 


and  //  is  the  thickness  of  the  inhomogeneous  layer,  p  the  den¬ 
sity,  X  and  p  the  Lame  constants.  For  topographic  irregularities 


<*(*)=  — M<^(x)i;(*  +y)>e*  >d*y  (4) 
(2nyJ 


where  f,  is  the  height  of  the  interface  with  respect  to  the  average. 

In  our  model  calculations  we  used  the  PREM  structure 
[Dziewonski  and  Anderson,  19*1]  to  compute  travel  time, 
geometric  spreading,  and  attenuation.  The  PKP  wave  has  a 
caustic  near  the  CMB  and  ray  theory  predicts  large  amplitudes 
in  its  vicinity.  We  compute  separately  the  contribution  from  both 
branches  AB  and  BC,  and  neglect  interference.  Near  the  caustic 
we  impose  an  amplitude  cutoff  equal  to  the  amplitude  one 
wavelength  away  from  the  cusp.  Changing  the  cutoff  amplitude 
to  the  amplitude  two  or  three  wavelengths  away  from  the  cusp 
did  not  appreciably  affect  the  result,  because  the  region  where 
caustics  develop  is  small  compared  to  the  total  scattering  region 
(several  tens  of  wavelengths).  The  observability  of  precursors  is 
not  due  to  developing  a  caustic  near  the  region  of  inhomo¬ 
geneities;  it  is  instead  an  effect  of  having  a  large  scattering 
region  at  the  CMB  and  of  being  the  first  arrivals  at  the  station. 

From  studies  of  large  seismic  arrays  [A'ing  et  ai,  1974: 
Doornbos,  1976|  it  was  possible  to  determine  for  some  cases  the 
location  of  scattering,  distinguishing  scattering  at  entrance  and 
exit  from  the  core.  For  our  data  this  is  not  possible  and  in  our 
calculations  we  include  both  contributions,  and  neglect  interfer¬ 
ence.  If  there  are  regions  for  which  only  scattering  at  entrance 
or  exit  contributes,  we  would  underestimate  the  strength  of  such 
inhomogenetties  by  roughly  one  half. 

For  volume  mhomogeneilies,  the  mean  square  amplitude  of 
the  scattered  waves  depend  on  three  independent  random  vari¬ 
ables,  p,  X,  and  p.  One  could  choose  another  set  of  independent 
random  variables  with  more  direct  physical  meaning,  such  as  p, 
p,  and  K  (compressibility);  p.  a  (P  wave  velocity),  and  (5  (5 
wave  velocity);  or  an  easier  seismological  interpretation  when 
X  =  p  by  a  and  pa  [l Vu  and  Aki,  1985a].  In  general,  it  is  not 
possible  to  determine  the  characteristics  of  each  variable 
independently.  Even  if  <t>(K )  is  known,  it  could  be  fitted  by  dif¬ 
ferent  combinations  of  strengths  and  correlations  among  these 
variables.  However,  if  we  consider  P  wave  scattering  in  the  for¬ 
ward  direction,  only  one  variable  is  involved:  a.  Backscattered 
P  waves  depend  on  the  P  wave  impedance,  Zp  =  pa.  P  wave 
scattering  at  90°  angle  depends  on  X  only.  Choosing  these  three 
(a,  Zp,  and  X)  as  independent  variables,  (3)  becomes 


5P  - ( cos©  +  cos2©)  +  — — 

a  Zp 

(cos©  -  cos2©)  -  — — — sin2©  (5) 

X  -i-  2p 


For  precursors  to  PKIKP  the  scattering  angle  varies  roughly 
between  10°  and  50°;  therefore  we  can  neglect  the  contribution 
from  fluctuations  of  P  wave  impedance.  Furthermore,  if  the 
inhomogeneities  are  due  to  temperature,  or  compositional  differ¬ 
ences,  leaving  Poisson's  ratio  constant  at  =  0.25,  such  that  we 
can  approximate  p  =  X,  6p  =  6X,  we  have 


50  =  -2-^-(cos2-y-  -  ysin2©)  (6) 

This  is  the  relation  we  use  to  interpret  volume  inhomogeneities 
near  the  CMB.  In  the  case  of  topographic  irregularities,  there  is 
only  one  random  variable  involved:  the  height.  From  (1),  (2), 
and  (4)  we  have 

<t>v(/f)  =  n<t>,(/ir)  (7) 

and  considering  a  thickness  of  200  km  for  the  inhomogeneous 

medium  a  correspondence  between  the  variation  of  a  and  the 
height  ^  is  given  approximately  by 

;s56—  (km)  (8) 

a 


The  average  scattering  angle  increases  with  decreasing  epi- 
ccntra!  distance.  For  a  given  epicentra!  distance  the  scattering 
angle  varies  by  no  more  than  about  10°.  Typically,  at  a  distance 
of  140°  the  average  scattering  angle  is  about  10°.  and  at  125°  the 
average  scattering  angle  is  about  50°.  If  the  spectral  density  of 
inhomogeneities  corresponding  to  the  scattering  angle  does  not 
change  significantly  for  variations  of  scattering  angle  of  10°, 
then  <t>(K)  can  be  assumed  constant.  Taking  the  ratio  of  ampli¬ 
tudes  of  precursor  /  PKIKP ,  from  ( 1 )  we  obtain 


<HK)  = 


<  lu\'l>pr“ 


d^ptapfleMp  )2  ( 
obs  _  '  calc 

nk,jA1R2ildS 


where  U  I ok,  is  the  observed  amplitude  of  PKIKP.  and  %pi‘kp 
and  Rpi,tp  are  the  geometric  spreading  and  transmission 
coefficients  relevant  to  PKIKP  calculated  using  ray  theory.  The 
effect  of  different  proposed  Q's  along  the  different  paths 
through  the  mantle,  outer,  and  inner  core  is  negligible  in  this 
case.  Fig.  6  shows  the  result  of  <t>(/C )  for  both  (1)  volumetric 
inhomogeneities  and  (2)  topograph- c  irregularities,  for  all  sta¬ 
tions.  Comparing  the  results  for  different  stations  we  find  no  evi¬ 
dence  for  regional  variations  of  the  strength  of  inhomogeneities 
of  greater  than  a  factor  of  three.  The  spectrum  can  be  approxi¬ 
mated  as  a  power  law  by  <t >(K)  =  I ‘'’(km4),  with  K  in 
units  of  km ~ 1  and  <t>0  =  (3.4  ±  1.5)  x  10^,  (4.3  ±  2.2)  x  10"7,  and 
p  =5.3  ±0.2,  6.8  ±0.3,  for  volumetric  inhomogeneities  and 
topographic  irregularities,  respectively.  This  is  different  from 
the  Gaussian  correlation  used  by  Haddon  and  Cleary  [  1974)  and 
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Fig.  6.  Spectra  <t>  (km4]  deduced  from  the  data  of  Figure  5  under 
the  assumption  of  (left)  D"  inhomogeneities  and  (right)  CMB 
topography.  The  line  represents  the  best  fit  power  law  which  has 
index  p  =  5.3  (left)  and  6.8  (right;.  The  spatial  wave  number 
I  XT  (km  ‘)  refers  to  the  inhomogeneities  (assumed  isotropic)  or 
to  the  topographic  features. 
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Doornbos  [1976],  Recently  F latte  and  Wu  11988]  inverted  for 
the  small-scale  inhomogeneities  beneath  NORSAR,  in  Norway, 
and  found  that  the  heterogeneities  below  the  lithosphere  are  best 
represented  by  a  power  law  spectrum  of  index  between  4  and  5. 

4.  Geophysical  Implications 

The  total  variance  of  the  P  wave  velocity  through  D"  and  the 
height  of  the  topography  at  the  CMB  are  obtained  from  the 
integration  of  <t>'J(K )  over  the  observed  domain  of  wavelengths 
(10-70  km).  The  result  for  topographic  irregularities  is 
280  ±  100-m  rms  height.  This  result  is  consisted  sith  the  esti¬ 
mate  by  Doornbos  [1978]  and  the  upper  limit  r  om  Me  nice 
[1986|  of  topographic  irregularities  of  a  few  hundred  meters 
height.  The  result  for  three-dimensional  heterogeneities  is 
0.5  ±0.1%  for  P  wave  velocity  variations.  It  should  be 
emphasized  that  our  result  determines  the  product  of  8a 2H ,  and 
we  are  assuming  that  H  =  200  km.  Relating  the  variations  of 
seismic  velocities  with  temperature  as  pSo’  =  a  ST .  with 
a  =  2  x  10* erg  cm  3o/(~ 1  (from  Jones  [1977]),  where 
9*  =  a2-— p2,  and  neglecting  variations  in  pressure)  we  obtain 

ST  =  700°K.  It  is  possible  then  that  these  heterogeneities  have  a 
thermal  origin. 

If  the  spectrum  of  the  correlation  function  of  topographic  irre¬ 
gularities  is  valid  for  larger  wavelengths,  then  one  can  predict 
the  increase  in  mean  square  height  expected  from  increasing  the 
domain  of  integration.  For  instance,  if  we  take  scales  up  to  1000 
km,  the  height  is  23  km,  which  is  somewhat  more  than  Creager 
and  Jordan  [1986]  suggest  (5  km)  for  the  height  of  the  irregular¬ 
ities  from  analysis  of  PKP  travel  times.  This  suggests  that  the 
true  spectrum  lies  below  a  simple  extrapolation  of  our  spectrum 
to  wavelengths  significantly  above  70  km.  It  may  be  worthwhile 
to  point  out  that  mechanisms  with  characteristic  scales  of  a  few 
hundred  kilometers  could  explain  this  behavior.  For  example,  a 
statistical  process  with  a  characteristic  scale  and  an  exponential 
correlation  function  has  a  flat  spectrum  at  long  wavelengths  and 
a  K  4  spectrum  at  short  wavelengths.  Examples  of  such  mechan¬ 
isms  could  be  convection  cells,  or  the  catchment  area  for  plume 
material.  Similar  extrapolations  for  volume  inhomogeneities  are 
not  possible  without  including  anisotropy.  It  seems  natural  to 
consider  the  thickness  of  the  D"  layer  as  the  upper  limit  for  the 
scales  where  the  inhomogeneities  can  be  considered  isotropic. 
This  is  approximately  3  times  the  upper  limit  of  our  analysis. 
Assuming  that  the  spectrum  of  the  correlation  function  of 
volumetric  inhomogeneities  is  valid  up  to  200  km,  a  value  of 
2.5%  for  the  P  wave  velocity  variation  is  predicted.  If  this  P 
wave  variation  has  a  thermal  origin,  variations  in  temperature 
would  have  to  be  ST  *=  3000°K.  New  results  from  laboratory 
experiments  have  indicated  that  the  temperature  change  across 
D"  is  about  2000°K  I Knittle  el  al.  1987],  so  that  a  thermal  origin 
for  our  observed  P  wave  variations  is  possible.  On  the  other 
hand,  the  derivative  of  P  wave  velocity  with  respect  to  tempera¬ 
ture  may  be  much  less  at  higher  pressure  [Ahrens  and  Hager, 
1987],  which  would  require  much  larger  temperature  variations 
to  explain  our  observations.  »*e  conclude  that  interpreting  our 
results  in  terms  of  thermal  variations  alone  would  have  to  await 
the  resolution  of  present  uncertainties  on  the  temperature  change 
across  D"  and  the  relation  of  temperature  change  to  P  wave 
velocity  change  al  high  pressure. 

It  is  also  possible  to  associate  our  observed  P  wave  velocity 
variations  with  chemical  heterogeneity.  This  heterogeneity  could 
arise  from  mixing  between  the  iron  core  and  mantle  silicates, 
resulting  in  formation  of  mantle  minerals  such  as  FcO  [ Knillle  el 


al.  1987],  Since  the  P  wave  velocity  of  FeO  alloys  is  expected 
to  be  substantially  smaller  than  that  of  the  mantle,  this  model 
would  require  very  small  inclusions  contributing  a  variable  but 
small  mass  fraction  to  a  200-km-thick  D".  Other  mechanisms 
include  subducted  crustal  slabs  [Ringwood,  1979.  Olson  el  al 
1987  a,  i>]  and  primary  chemical  layering  [Jordari  and  Creager, 
1987], 


5.  Conclusions 

By  analyzing  the  characteristics  of  precursors  for  all  available 
source -receiver  paths  from  GDSN  data,  we  have  shown  that  the 
spectrum  of  inhomogeneities  satisfies  a  power  law  distribution 
of  index  5.3  ±  0.2  for  volumetric  inhomogeneities  and  6  8  ±  0.3 
for  topographic  irregularities.  These  inhomogeneities 
correspond  to  variations  of  280  ±  100-m  rms  height  if  scattering 
is  due  only  to  topography,  or  approximately  0.5  ±  0.1%  in  P 
wave  velocity  if  it  is  due  only  to  volumetric  inhomogeneities.  At 
present,  we  can  not  distinguish  between  these  two  types  of  mho- 
mogeneities.  We  find  evidence  for  similarities  of  the  properties 
of  inhomogeneities  (P  wave  speed  variation  in  D"  or  rms  height 
of  CMB  topography)  between  several  different  geographical 
regions  within  a  factor  of  three  in  the  strength  of  the  inhomo¬ 
geneity.  It  is  possible  at  present  to  consider  these  variations  as 
arising  from  either  thermal  or  chemical  heterogeneities,  but  this 
information,  in  combination  with  further  constraints  on  material 
properties  from  laboratory  experiments  and  progress  in  geo¬ 
dynamic  simulations,  will  significantly  limit  allowable  models 
of  the  core-mantle  boundary  region. 

Appendix 

Here  we  present  briefly  the  theory  of  wave  scattering  in  ran¬ 
dom  media  used  in  this  study.  Two  cases  are  analyzed:  scatter¬ 
ing  by  volumetric  inhomogeneities  and  topographic  irregulari¬ 
ties. 


Volume  Inhomogeneities 


We  use  a  theory  developed  by  Chernov  [1960]  for  acoustic 
waves  and  extended  to  elastic  waves  by  Knopoff  and  Hudson 
[1964]  and  Wu  and  Aki  [1985a,  b].  This  method  was  also  used 
for  studies  of  precursors  to  PK1KP  by  Haddon  and  Cleary 
[1974]  and  Doornbos  [1976], 

The  equation  of  motion  for  a  general  inhomogeneous  isotro¬ 
pic  elastic  medium  is 

p-^j«(*,0-v|\v  ■«(*,!)]  - 


■M 


Vu(s,  t )  +  Vu(s 


•'>1 


t  =  /(*.') 


(Al) 


or 


Lu  =  /  (A2) 

where  p  is  the  density,  X  and  p  are  the  Lame  constants,  u  is  the 
displacement,  /  is  the  body  force,  and  L  a  linear  operator  which 
depends  on  the  parameters  of  the  medium.  For  a  radially  sym¬ 
metric  Earth  model,  p0,  Xo  and  p<>  vary  smoothly  with  radius, 
and  (A2)  can  be  written  as 

Loum  =  f  (A3) 

with  solutions  given  by 

u(0'0t,/)  =  jGU.x'.r)  */(*'.  OdV  (A4) 

where  the  Green’s  function  C(x,x',l)  is  obtained  from  ray 
theory  [ Aki  and  Richards.  1980],  and  *  means  convolution  ui 
lime.  In  the  presence  of  inhomogeneilies  such  that 
p  =  po  +  5p,  X  =  X<|  +  5X,  and  p  =  Pe  +  op,  the  equation  of 
motion  assuming  no  body  forces  becomes 
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IS.Ofri 


L  cM  =  Qu 


where 

Q  =  V(5XV  )  +  V  j^5n(V  +  Vr)j  -  5p-|^ 


(A5) 

(A6) 


depends  only  on  the  inhomogeneities. 

Decomposing  u  into  the  incident  wave,  «(0),  plus  the  scattered 
wave,  un\  and  assuming  that  l«(1)l  «  lu(0)l  (Bom  approxima¬ 
tion),  the  equation  for  the  scattered  wave  is 

Louw  --  Qum  (A7) 

which  has  a  solution 


«<'>(*,  r)  =  jGOr.jc'.r)  *  Qu'a\x',t)d2x'  (A8) 

To  calculate  the  average  intensity  of  the  scattered  P  wave 
over  all  possible  realizations  of  the  random  medium,  we  take  the 
ensemble  average  of  the  square  amplitude  and,  after  some  alge¬ 
bra,  find 


<  lu(1)l2>  =  jkijcMBA2(x’)R2tt>(K)d2x'  (A9) 

where  k  is  the  wave  number  of  the  wave  field,  A  is  the  product 
of  the  geometric  spreading  coefficient  from  source  to  the  scatter¬ 
ing  point  and  the  geometric  spreading  coefficient  from  the 
scattering  point  to  the  receiver,  R  contains  all  the 
transmission/reflection  coefficients  and  attenuation  factors, 
K  =  k(n  -  ri’),  and  n  and  ri'  are  the  direction  of  the  incident  and 
scattered  wave  at  the  scattering  point.  Assuming  0  to  be  the 
scattering  angle,  and  //  the  thickness  of  the  inhomogeneous 
layer  we  have 

<t>(K )  =  — ~rf  <  &3(x)5C)(x  +  y)  >  ex  yd2y  (AlOa) 
(2  ny 

5d(x)  =  -^-cos©  -  -t—— - — cos2©  (AlOb) 

p  X  +  2p  X  +  2p 

From  the  analysed  earthquakes  which  have  magnitudes  5.7 
<  mb  <  6.7,  the  duration  of  the  PKIKP  wave  is  between  4 
and  8  s,  and  the  dominant  frequency  is  around  1  Hz.  To  model 
these  charactcrisitics  we  assumed  a  source  of  monochromatic 
waves  with  a  frequency  of  1  Hz  and  a  Gaussian  envelope  of  4  s 
duration.  Changing  the  source  duration  from  4  to  8  s  does  not 
affect  the  result  appreciably. 

Topographic  Irregularities 

No  exact  theory  of  wave  scattering  for  general  rough  surfaces 
is  available.  However,  effective  approximate  methods  have 
been  developed  for  some  cases.  The  small  perturbation  method 
which  is  valid  for  roughness  with  small  heights  and  slopes  was 
applied  to  elastic  media  by  Kennelt  [1972]  and  used  to  study 
precursors  to  PKIKP  by  Doornbos  [1978).  On  the  other  hand, 
the  Kirchhoff-Helmholtz  method  is  valid  for  large  heights  with 
smooth  slopes  and  was  used  in  seismology  by  Scott  and  Helm- 
berger  [1983],  who  studied  body  wave  reflections  from  moun¬ 
tains  and  spall  from  nuclear  blasts;  Frazer  and  Sen  [1985]  and 
Sen  and  Frazer  [1985]  used  it  to  study  reflections  in  a  laterally 
inhomogeneous  multilayered  medium.  Here  we  present  a  brief 
derivation  of  the  Kirchhoff-Helmholtz  method. 

Let  V  be  a  volume  and  dV  its  boundary  with  outward  point¬ 
ing  unit  normal  rri.  Let  /,  and  f  i  be  some  distribution  of  force 
densities,  and  let  u , ,  X| ,  iz  2,  and  t2  be  the  displacements  and 
stresses  due  to  /,  and  /2  respectively.  The  equation  of  motion 
in  the  frequency  domain  is 


-  pa>2u,  =  V  t,  +  /,  (i  =  1,2)  (All) 


Integrating  (11)  with  i  =  1  over  a  volume  V  after  dot  product 
with  u  2.  permuting  indices  and  subtracting  from  each  other,  and 
using  the  divergence  theorem,  we  obtain 


|(/2  “i -ft  u7)dV  =  |  rri  ■  (r,  Uj-t2'  u{)dS  (A12) 


Here  f  \  and  / 2  are  in  general  arbitrary;  by  choosing  /,  to  be 


nonzero  outside  of  V  and  dV,  and  /2  =  dSfr  -  *2)  where 
x-i  e  V ,  we  obtain 

d“,(*i)=  f  m -(ti  •  u2-t2- u,)dS  (A13) 

oV 

Considering  a  medium  with  a  material  discontinuity,  the 
boundary  dK  is  taken  as  a  surface  infinitely  close  to  the  material 
discontinuity  from  the  receiver  side.  Using  ray  theory  to  calcu¬ 
late  Ui,  t|,  u2.  and  t2,  the  expression  for  the  displacement,  drop¬ 
ping  the  subscript  1 ,  becomes 


where  I  is  the  illuminated  surface  of  the  CMB;  A  and  R  are  the 
same  as  in  (A9),  T  is  the  travel  time,  and  y  depends  on  the  slope 
of  the  interface  and  physical  parameters  of  both  media  as  fol¬ 
lows; 


P'S' 


(A15) 


Cp  =  (X  -  2pri  ■  ri'Xn'  -  ri). 


C,  =  (-jpj"  -  2pri')s  ■  n  +  r(yfLf"  ’  -  *•)  (A16) 

where  ri  is  the  direction  of  the  incident  wave,  ri',  s',  P'P'  and 
P’S'  are  the  directions  of  displacement  and  transmission 
coefficients  of  transmitted  P  and  S  waves,  respectively,  and  s' 
is  the  direction  of  propagation  of  the  transmitted  S  wave  satisfy¬ 
ing  s'  •  s  =  0. 

Let  C,(r )  be  the  height  of  the  interface  with  respect  to  the  aver¬ 
age.  For  small  heights,  A  is  considered  independent  of  and 
keeping  only  the  first  term  of  the  expansion  of  R  and  y  in 
powers  of  V£,  and  T  in  powers  of  the  ensemble  average  of 
the  square  amplitude  takes  the  form 


<  I  u  l 2  > 


a26(K)  +  k2iW(K  )!■ 


d2r  (A17) 


where  ft  depends  only  on  the  geometry  and  elastic  constants  of 
both  media. 


Y(/0=  ’  f<^(x)^  +y)>e*  ’d2 y 

(271) 


(AI8) 

(A19) 


where  the  prime  means  derivative  with  respect  to  C,,  and  sub¬ 
script  0  refers  to  ^  =  0. 

The  first  term  of  the  right  hand  side  of  (A  17)  contributes  to 
the  specular  direction,  because  K  =  0  is  just  Snell’s  law.  The 
second  term  describes  finite  angle  scattering  waves  and  has  a 
similar  form  to  that  for  the  volumetric  inhomogeneities  obtained 
in  (A9),  except  for  the  geometry  factor  ft.  The  validity  of  this 
method  is  limited  to  topographic  irregularities  with  radius  of 
curvature  larger  than  the  wavelength,  so  that  locally  plane 
reflection/transmission  coefficients  can  be  used,  or  h  «  L  and 
X  «  L2/8ri,  where  h  and  L  are  the  characteristic  height  and 
length  of  the  irregularity  respectively,  and  X  is  the  seismic 
wavelength.  In  the  application  of  this  method  to  PKP  waves 
where  scattering  occurs  at  the  CMB,  the  wavelength  of  P  waves 
of  1  Hz  is  typically  about  13  km,  so  that  this  method  can  be  used 
for  irregularities  with  heights  of  up  to  1  km  and  wavelengths  not 
smaller  than  10  km. 
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Small-Scale  Structure  in  the  Lithosphere  and  Asthenosphere 
Deduced  From  Arrival  Time  and  Amplitude  Fluctuations 

at  NORSAR 


Stanley  M.  Flatte  and  Ru-Shan  Wu 

Physics  Department  and  Institute  of  Tectonics,  University  of  California,  Santa  Cruz 


We  analyze  the  pattern  of  phase  and  amplitude  variations  of  seismic  waves  across  the  NORSAR  array  on  a 
statistical  basis  in  order  to  determine  the  statistical  distribution  of  heterogeneities  under  NORSAR.  Impor¬ 
tant  observables  that  have  been  analyzed  in  the  past  are  the  phase  (or  travel  time)  and  log  amplitude  vari¬ 
ances  and  the  transverse  coherence  functions  (TCFs)  of  phase  and  amplitude  fluctuations.  We  propose  and 
develop  the  theory  and  methods  of  using  other  observables  to  reduce  tire  degree  of  nouniqueness  and 
increase  the  spatial  resolution  of  the  analysis.  Most  important  are  the  angular  coherence  functions  (ACFs), 
which  characterize  quantitatively  the  change  in  the  pattern  of  fluctuations  across  the  array  from  one  incom¬ 
ing  angle  (or  beam)  to  another  and  which  have  a  different  sensitivity  to  the  depth  distribution  of  hetero¬ 
geneities  than  the  TCFs.  A  combination  of  the  ACFs  and  TCFs  allows  estimation  of  the  power  spectra  of 
the  P  wave  speed  variations  under  the  array  as  a  function  of  depth.  We  use  data  for  phase  fluctuations 
from  104  incident  beams  and  amplitude  fluctuations  from  18S  beams  with  2-Hz  center  frequency  at  NOR¬ 
SAR  to  calculate  the  three  ACFs  and  three  TCFs  (of  phase,  log  amplitude,  and  their  cross  coherence).  The 
measured  ms  travel  time  fluctuation  is  0.135  s,  and  the  rms  log  amplitude  fluctuation  is  0.41.  The  half¬ 
coherence  widths  of  the  ACFs  are  3°  for  log  amplitude  and  9°  for  phase.  The  half-coherence  widths  of  tile 
TCFs  ate  18  km  for  phase  and  less  than  the  minimum  separation  between  the  elements  of  the  array  for  log 
amplitude  In  order  to  account  for  these  features  of  the  data,  we  adopt  a  two- overlapping -layer  model  for 
lithospheric  and  asthenospheric  heterogeneities  underneath  NORSAR,  with  spectra  drat  are  band-limited 
between  the  wavelengths  of  5.5  and  1 10  km.  Our  best  model  has  an  upper  layer  with  a  flat  power  spec¬ 
trum  extending  from  the  surface  to  about  200  km,  and  a  lower  layer  with  a  K~*  power  spectrum  extending 
from  15  to  250  km.  The  latter  spectrum  corresponds  to  an  exponential  correlation  function  with  scale 
larger  than  the  observation  aperture  (110  km).  The  rms  P  wave  speed  variations  lie  in  the  range  1-4%. 
The  small  scale  heterogeneities  may  be  attributed  to  clustered  cracks  or  intrusions;  the  larger-scak  wave- 
speed  heterogeneities  are  temperature  or  compositional  heterogeneities  that  may  be  related  to  chemical  dif¬ 
ferentiation,  or  dynamical  processes  in  the  boundary  layer  of  mantle  convection. 


1 .  Introduction 

Observations  of  direct  P  wave  amplitude  and  arrival  lime 
fluctuations  [AJa,  1973;  Capon,  1974;  Capon  and  Berteussen, 
1974;  Berteussen  et  at.,  1975 a,  b\  Powell  and  Meltzer,  1984], 
coda  strength  [Aki,  1981;  Sato,  1982;  Wu  and  Aki,  1985a,  b], 
and  attenuation  by  scattering  [ Aki  and  Chouet,  1975;  Aki,  1980; 
Wu,  1982;  Sato,  1982]  have  all  been  used  in  attempts  to  deter¬ 
mine  some  statistical  characteristics  of  small-scale  structure  in 
the  Earth.  Statistical  analyses  of  amplitude  and  arri'"l  time 
fluctuations  have  previously  involved  the  variances  of  log  ampli¬ 
tude  and  arrival  time,  the  covariance  of  log  amplitude  and 
arrival  time,  and  the  coherence  functions  of  log  amplitude  and 
arrival  time  as  functions  of  spatial  separation  along  the  Earth's 
surface  (the  transverse  coherence  functions,  or  TCFs)  [Ah', 
1973;  Capot  1974;  Capon  and  Berteussen,  1974;  Berteussen  et 
ai,  1975 a,  b\ .  ’Hand Meltzer,  1984], 

Theoretical  analyses  of  the  observations  assume  a  particular 
medium  structure  and,  by  various  techniques,  compare  theoreti¬ 
cal  predictions  of  seismic  wave  properties  to  the  observations. 
Previous  seismic  wave  analyses  followed  Chernov  [1960],  who 
used  the  Rytov  and  Fresnel  approximations  to  connect  wave 
fluctuations  to  medium  variations.  His  expressions  involving 
general  random  media  in  the  space  domain  are  formidable,  but 
he  evaluated  them  explicitly  for  a  statistically  uniform  and  iso- 
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tropic  random  medium  with  a  Gaussian  correlation  function;  hit 
resulting  expressions  have  been  used  by  previous  workers  to 
obtain  strengths  and  scale  sizes  of  the  P  wave  speed  variations 
in  the  lithosphere  by  modeling  it  as  a  uniform  layer  of  a  single- 
scale  random  medium. 

This  paper  has  several  purposes;  first,  to  point  out  that  the  log 
amplitude  and  arrival  time  fluctuations  at  large  seismic  arrays 
like  NORSAR  and  LASA  have  further  important  statistical 
information  that  has  not  been  utilized  in  previous  analyses,  espe¬ 
cially  the  coherence  function  of  waves  from  different  sources, 
and  hence  different  incoming  directions;  second,  to  present  a 
statistical  analysis  of  NORSAR  data  that  includes  this  new 
information;  third,  to  present  results  from  the  modem  theory  of 
wave  propagation  through  random  media  (WPRM)  that  are 
based  on  the  parabolic  wave  equation  approach  and  are  some¬ 
what  simpler  than  Chernov’s  because  they  are  formulated  par¬ 
tially  in  the  spectral  domain  and  can  accommodate  easily  any 
model  of  the  medium  spectrum;  and  fourth,  to  present  a  model 
of  the  inhomogeneities  in  the  lithosphere  and  asthenosphere 
under  the  NORSAR  array  that  is  consistent  with  the  available 
data. 

Before  discussing  the  most  interesting  of  the  new  information 
available,  it  is  important  to  point  out  that  there  are  three  vari¬ 
ances  involved  in  arrival  time  and  log  amplitude  fluctuations  at  a 
given  seismic  frequency;  the  variances  of  arrival  time  and  log 
amplitude,  and  the  covariance  between  the  two.  Each  of  these 
variances  leads  to  a  coherence  function  of  any  variable  being 
investigated;  for  example,  previous  analyses  [Aki,  1973;  Capon, 
1974;  Capon  and  Berteussen,  1974;  Berteussen  et  al.,  1975 a,  b; 
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Powell  and  Meltzer,  1984]  have  used  the  arrival  time  and  log 
amplitude  transverse  coherence  functions,  but  we  will  be  analyz¬ 
ing  here  for  die  first  time  the  transverse  cross-coherence  func¬ 
tion  between  arrival  time  and  log  amplitude. 

The  main  new  information  that  can  be  obtained  from  an  array 
that  detects  waves  arriving  from  many  directions  ("beams")  is 
the  set  of  three  angular  coherence  functions  (ACFs):  that  is,  the 
coherence  functions  of  arrival  time,  log  amplitude,  and  their 
cross  coherence,  all  as  a  function  of  the  angle  between  two 
incoming  directions.  Since  the  difference  in  direction  between 
two  beams  may  be  as  small  as  1°,  the  information  at  a  given 
receiver  probes  (in  a  statistical  sense)  inhomogeneities  that  are 
quite  small:  of  the  order  of  1  km  at  60- km  depth.  The 
transverse  coherence  functions  from  a  coarse  array  with  -10-km 
spacing  cannot  probe  to  scales  smaller  than  10  km. 

Data  on  arrival  time  and  amplitude  fluctuations  of  the  first- 
arrival  teleseismic  P  wave  signal  are  analyzed  in  a  different 
fashion  than  data  involving  large-angle  scattering  of  wave 
energy.  The  difference  involves  the  realization  that  the  full 
wave  equation  need  not  be  solved,  but  rather  a  simpler  equation: 
the  parabolic  wave  equation  (PWE)  that  adequately  treats  waves 
in  a  narrow  angular  cone.  (Note  that  if  coda  is  interpreted  as 
large-angle  scattered  waves,  it  cannot  be  treated  in  tlie  same 
way.)  One  consequence  of  using  the  PWE  is  that  the  theoretical 
formulas  are  easily  expressed  in  terms  of  an  integral  along  the 
unperturbed  ray,  rather  than  as  a  volume  integral  over  all  of 
space.  Section  2  briefly  introduces  the  modem  theory  of  WPRM 
theory  based  on  the  PWE  and  "weak"  fluctuations.  Section  3 
specifically  describes  the  theory  of  the  angular  coherence  func¬ 
tions. 

The  NORSAR  data  used  in  our  analysis  consist  of  the  travel 
time  anomalies  [Berteussen,  1974]  at  22  subarrays  for  104 
beams;  and  the  log  amplitude  fluctuations  [Berteussen  and 
Husebye ,  1974]  at  22  subarrays  for  185  beams.  All  the  beams 
have  incoming  directions  within  a  35°  cone  around  the  zenith. 
The  smallest  distance  between  subarrays  is  about  10  km,  and  the 
largest  is  about  1 10  km.  The  data  used  were  filtered  for  1-3  Hz, 
so  the  nominal  frequency  is  2  Hz.  It  is  important  to  realize  that 
the  data  that  we  analyze  are  insensitive  to  inhomogeneities  with 
wavelengths  greater  than  about  100  km  because  of  the  finite  size 
of  the  array  and  are  also  insensitive  to  inhomogeneity 
wavelengths  less  than  about  5  km  because  of  two  effects:  first 
the  data  are  averages  over  subarrays  that  are  7  km  in  diameter, 
and  second,  the  wavelengths  of  the  seismic  waves  are  about  4 
km.  Section  4  describes  the  data  analysis  and  presents  results 
for  coherence  functions  at  NORSAR. 

It  does  not  take  a  sophisticated  theory  to  draw  some  dramatic 
conclusions  from  this  new  analysis  of  NORSAR  data.  Briefly, 
the  log  amplitude  ACF  drops  rapidly,  reaching  a  value  of  0.5  at 
an  angle  of  2°,  followed  by  a  more  gradual  drop  to  0.1  at  10°.  In 
contrast,  the  arrival-time  ACF  drops  much  more  slowly,  imply¬ 
ing  much  larger -scale  structure.  Yet  all  the  TCFs  at  first  glance 
have  scales  in  the  range  of  10-20  km.  Section  5  contains  the 
quantitative  comparison  of  weak  scattering  models  with  the 
NORSAR  data.  We  have  found  that  medium  spectral  models 
that  are  homogeneous  in  depth  (that  is,  the  spectrum  does  not 
change  with  depth)  down  to  a  cutoff  depth  between  0  and  500 
km  cannot  fit  all  the  data  at  once;  in  particular,  such  models, 
which  can  fit  the  TCFs  by  themselves,  cannot  simultaneously  fit 
the  shapes  of  the  ACFs.  We  suggest  a  two-overlapping-layer 
model  in  which  the  medium  spectrum  is  characterized  by  a 
power  law  K~*  over  the  sensitive  wave  number  band 


(5<  X  clOOkm).  An  upper  layer  with  significant  small-scale 
structure  (p  =  0)  spans  0<  z  <  200km  and  a  lower  layer  with 
strength  concentrated  at  large  scales  (p  =  4)  spans 
15  <z  <  250  km.  The  P  wave  rms  velocity  variations  in  our 
model  vary  between  1  and  4%.  A  specific  prediction  of  the 
model  that  we  suggest  is  a  rapid  drop  in  the  log  amplitude  TCF 
for  separations  of  a  few  kilometers,  followed  by  a  more  gradual 
drop  with  a  scale  of  tens  of  kilometers.  -  -  have  checked  this 
prediction  by  determining  the  TCF  of  log  amplitude  from  the 
individual  stations  at  NORSAR,  using  13  nuclear  explosion 
events.  We  find  consistency  with  our  prediction,  giving  us  some 
confidence  in  our  interpretation. 

2.  Theory  Of  WPRM 

Many  aspects  of  seismic  wave  propagation  in  the  Earth  involve 
examples  of  a  more  general  branch  of  science;  namely,  wave 
propagation  through  random  media  (WPRM).  Progress  in  this 
field  has  involved  coupling  of  the  theory  of  WPRM  with  investi¬ 
gation  into  the  detailed  character  of  the  (generally  fluid)  random 
medium.  At  least  two  theoretical  approaches  may  be  taken:  first, 
the  problem  may  be  approached  from  a  deterministic  point  of 
view  in  which  the  analysis  is  carried  out  for  particular  complex 
media.  Numerical  simulation  is  the  extreme  of  this  point  of 
view.  Second,  the  problem  may  be  treated  statistically  from  the 
outset.  In  this  approach,  one  assumes  a  spectral  model  for  the 
medium  and  attempts  by  analytical  means  to  predict  the  statisti¬ 
cal  behavior  of  the  propagating  wave  field.  Here  we  take  this 
second  point  of  view. 

The  modem  theory  of  WPRM  may  be  said  to  have  begun  in 
the  late  1940s  and  early  1950s  [Bergmann,  1946;  Mintzer,  1953; 
Chernov,  1960;  Tatarskii,  1971]  when  researchers  used  pertur¬ 
bation  techniques  to  develop  general  formulas  for  propagation 
through  weak  fluctuations.  These  formulas  involved  die  Bom 
approximation  and  hence  a  volume  integral  over  the  medium 
inhomogeneities.  In  comparing  with  experiment  they  assumed 
that  the  medium  was  characterized  by  a  Gaussian  correlation 
function,  but  unfortunately,  no  natural  medium  is  known  with 
this  property. 

The  problem  of  radio  wave  propagation  in  and  through  the 
Earth's  ionosphere  has  been  of  interest  since  the  1940s.  This 
field  made  a  crucial  contribution  to  WPRM  through  the  work  of 
Leontovich  and  Fock  [1946]  and  Fock  [1950],  who  introduced 
the  parabolic  equation  method,  which  treats  waves  that  are  con¬ 
centrated  within  a  small  angular  region  around  the  direction  of 
propagation.  Nearly  all  subsequent  analytic  work  has  used  the 
parabolic  equation  as  a  starting  point,  which  restricts  validity  to 
waves  with  directions  confined  within  a  cone  of  full  angle  about 
30°.  A  great  advantage  of  the  parabolic  equation  is  that  it  results 
in  formulas  requiring  line  integrals  along  deterministic  rays, 
rather  than  volume  integrals. 

A  major  step  forward  in  WPRM  was  taken  in  the  1960s,  in 
response  to  developing  understanding  of  Kolmogorov's  pioneer¬ 
ing  characterization  of  the  p  =  11/3  power  law  spectrum  of 
homogeneous,  isotropic  turbulence  (HIT).  A  theory  for  wave 
propagation  through  HIT  was  developed,  primarily  by  Soviet 
workers,  and  successfully  applied  in  the  case  of  weak  fluctua¬ 
tions  to  light  transmission  through  the  atmosphere,  making  use 
of  the  parabolic  equation.  Their  work  was  summarized  in  the 
influential  book  by  Tatarskii  [  1971]  and  a  review  article  by  Pro¬ 
khorov  et  al.  [1975],  and  more  recent  results  in  this  field  are 
covered  in  the  review  articles  by  Ishimaru  [1977],  Fante  [1980], 
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and  Tatarskii  and  Zavorotnyi  [1980],  as  well  as  in  the  book 
edited  by  Strohbehn  [  1978]. 

It  is  now  understood  that  weak  fluctuation  theory  applies 
when  the  variance  of  log  amplitude  is  small  compared  with 
unity.  We  realize  that  many  seismic  data  sets  have  large  inten¬ 
sity  variations,  but  in  this  paper  we  consider  only  those  for 
which  the  log  amplitude  variance  is  small.  For  example,  the 
NORSAR  data  set  of  subanay  amplitude  variations  has  a  log 
amplitude  variance  of  0.2.  We  defer  discussion  of  the  problem 
of  strong  fluctuations  to  the  future. 

Seismic  rays  often  turn  through  large  angles,  so  that  the  para¬ 
bolic  approximation  is  not  clearly  appropriate.  However,  if  the 
wave  fluctuations  are  caused  by  medium  variations  that  are 
confined  to  regions  near  the  source  and  receiver  (within  a  few 
hundred  kilometers),  then  the  rays  are  not  changing  their  angles 
significantly,  and  the  parabolic  approximation  is  valid. 

Work  on  fluctuations  in  sound  transmission  through  the  ocean 
during  the  1960s  was  largely  a  misguided  attempt  to  graft  the 
concept  of  HIT  onto  ocean  variability  [ Tatarskii ,  1971].  By  the 
mid-1970s,  oceanographers  had  identified  internal  waves  as  the 
most  important  source  of  variability  on  lime  scales  from  a  few 
minutes  to  a  day  [Garrett  and  Munk  1975], 

The  ocean  internal  wave  medium  provided  a  challenge  to 
those  interested  in  WPRM,  particularly  in  strong  fluctuations.  A 
significant  response  to  this  challenge  was  developed  over  the 
late  1970s.  The  first  success  in  this  area  was  achieved  for  weak 
fluctuations  by  Munk  and  Zachariasen  [1976],  whose  absolute 
calculations  of  variances  in  phase  and  log  amplitude  from  inter¬ 
nal  wave  effects  were  within  a  factor  of  2  of  the  available  exper¬ 
imental  results.  Strong  fluctuations  were  another  matter. 

A  further  basic  step  forward  took  place  between  1975  and 
1977,  with  the  consideration  of  the  arbitrary  field  of  wave  speed 
in  Feynman's  path  integral  as  a  statistical  random  medium 
[Dashen,  1979;  Flattd  el  al,  1979),  in  combination  with  the 
understanding  of  the  ocean  internal  wave  field  as  the  source  of 
medium  fluctuations.  These  ideas  led  to  successful  comparisons 
with  a  number  of  ocean  acoustics  experiments,  summarized  in 
the  monograph  by  Flatle'  et  al.  [1979],  the  review  article  by 
Flattd[  1983],  and  subsequent  publications  [Dashen  et  al.,  1985; 
Reynolds  et  al.,  1985;  Flattd  and  Stoughton,  1986;  Stoughton  et 
al,  1986;  Flattd.  et  al,  1987a.  b\. 

With  this  background  in  mind,  we  now  apply  parabolic  equa¬ 
tion,  weak-scattering  theory  to  P  wave  teleseismic  propagation. 
We  take  the  point  of  view  that  the  seismic  waves  arrive  a  few 
hundred  kilometers  below  an  array  like  NORSAR  in  essentially 
plane  wave  form:  that  any  effect  due  to  inhomogeneities  near  the 
source  is  too  large  scale  to  be  observed  by  a  100-km  array  and 
that  the  very  deep  mantle  has  imposed  no  significant  structure  on 
the  wave  that  has  scales  less  than  the  100-km  size  of  the  array. 
The  reason  that  the  inhomogeneities  near  the  source  have  little 
effect  is  that  only  a  very  small  range  of  initial  angles  at  the 
source  is  seen  by  the  array,  and  the  spreading  of  the  ray  tube 
magnifies  any  small-scale  fluctuations  of  the  wave  front  into 
large-scale  fluctuations  at  the  receiver. 

Furthermore  we  assume  that  the  effect  on  the  arrival  time  and 
amplitude  fluctuations  of  inhomogeneities  within  a  few  hundred 
kilometers  depth  below  the  array  results  mainly  from  wave 
scattering  within  a  reasonably  narrow  cone.  We  realize  that 
coda  can  represent  scattering  at  larger  angles;  however,  the  first 
few  seconds  of  P  wave  arrivals  are  restricted  to  small  angles. 

Since  NORSAR  data  (by  subanay)  show  a  log  amplitude 
variance  of  0.2,  we  use  weak  fluctuation  theory.  It  has  been 


Fig.  1.  Schematic  of  P  wave  propagation  through  the  upper  mantle, 
lithosphere,  and  crust  to  a  surface  seismic  array.  The  heavy  solid 
straight  line  is  the  incident  wave  front,  the  light  solid  horizontal  line  is 
the  wave  front  arriving  at  the  receiving  array  in  the  absence  of  fluctua¬ 
tions.  The  wavy  line  could  represent  either  the  wave  front  arrival  time 
accros  the  amy  or  the  intensity  at  each  point  in  the  arriving  wave  front 

established  that  even  if  intensity  fluctuations  are  large,  if  they 
are  averaged  over  a  small  source  (or  receiver)  region,  such  that 
the  residual  (larger  scale)  intensity  fluctuations  are  small,  then 
weak- fluctuation  theory  can  be  applied  [Wong  et  al,  1978]. 
There  are  probably  aspects  of  this  averaging  in  the  NORSAR 
data  because  of  the  7 -km  width  of  a  typical  subarray  and 
because  of  finite  wavelength  effects,  which  result  in  the 
observed  log  amplitude  variance  of  0.2. 

Finally,  we  neglect  the  small  deterministic  refraction  in  the 
few  hundred  kilometers  below  the  array;  in  the  absence  of  inho¬ 
mogeneities,  the  wave  fronts  remain  plane  in  the  same  direction. 
(See  Figure  1.) 

We  now  state  the  results  of  weak-scattering  theory  under  the 
above  restrictions.  First,  define  the  medium  P  wave  speed  as  a 
function  of  position: 

C(x)  -  C.[l+|i(x)]  (1) 

where  C„  is  a  reference  deterministic  speed  that  might  depend 
on  depth  z ,  and  p  ( x )  is  a  random  function  with  zero  mean  that 
represents  inhomogeneities.  The  statistics  of  p  are  assumed  to 
be  quasi-homogeneous.  That  is,  we  can  define  a  spectrum  at 
depth  z,  and  this  spectrum  is  allowed  to  change  slowly  with 
depth.  In  this  case,  the  character  of  p  (  x )  is  described  by  a 
three-dimensional  spectrum  W  (  K ,  z  ),  such  that 

<p(x)p(x')>  =  jd}KW(K,z)«“  (2) 

where  we  have  assumed  that  the  spectrum  depends  on  a  three- 
dimensional  wave  vector  K  and  depth  z ,  and  the  angle  brackets 
represent  an  average  over  the  statistical  ensemble  of  random 
media.  Since  the  spectrum  does  not  depend  on  transverse  posi- 
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tion  and  depends  only  slowly  on  depth  I ,  the  medium  correla¬ 
tion  function  <  p.  (  x )  |X  (x ')  >  depends  only  on  the  difference 
vector  x—  x,  except  for  a  slow  dependence  on  z  which  is 
implicit.  Note  the  normalization  of  the  spectrum  is  such  that 

<P*>  -  jd'KWiK.t)  (3) 

where  <  |i.2>  is  the  variance  of  the  fractional  change  in  wave 
speed  due  to  the  random  inhomogeneities;  it  can  depend  on 
depth  z . 

Consider  a  given  frequency  <o  in  the  incoming  wave.  That 
frequency  has  wave  number  k  =  a >  /C„ .  Let  the  complex  wave 
function  arriving  at  the  center  of  the  array  be 

4 >mAei"V'  (4) 

where  %  is  the  wave  function  that  would  arrive  in  the  absence 
of  die  random  inhomogeneities  and  where  the  amplitude  A  and 
phase  0  are  defined  in  this  relation.  Define  the  log  amplitude  as 
u,  so  that 


bi(4'/*K<))  =  u  +  i  0  (5) 

It  is  known  in  the  case  of  weak  scattering  confined  to  an  narrow 
cane,  in  an  intrinsically  dispersion- free  medium,  that  the  arrival 
time  fluctuation  t  and  the  phase  fluctuation  6  are  in  one-to-one 
correspondence: 

$  =  on 

We  now  point  out  that  for  a  specific  realization  of  inhomo¬ 
geneities,  u  and  $  are  functions  of  both  the  receiver  position  x 
and  the  source  position,  represented  by  a  unit  vector  0  describ¬ 
ing  the  direction  of  the  incoming  plane  wave.  (See  Figure  1 .) 

Note  first  that  we  have  made  the  approximation  that  the  ran¬ 
dom  medium  is  "horizontally  homogeneous";  that  is,  the  spec¬ 
trum  is  independent  of  horizontal  position.  We  also  consider 
waves  arriving  within  30®  of  the  zenith,  so  we  avoid  some  sim¬ 
ple  geometrical  corrections  that  are  important  only  for  highly 
slanted  rays.  We  can  split  the  three-dimensional  wave  number 
K  into  a  component  along  the  z  direction  K,  and  a  two- 
dimensional  vector  Kt  transverse  to  z;  that  is,  Kr  is  in  the  x  -y 
plane.  Then  the  medium  spectrum  can  be  expressed  as 

W(K,z)*W(.K,,Kt,z)  (6) 

With  the  above  definitions  and  approximations,  the  result 
[Afunk  and  Zachariasen.  1976;  Flattd  el  ai.%  1979;  S.M.  FI  and 
and  T.  Moody,  manuscript  in  preparation,  1987]  for  the  vari¬ 
ances  of  the  seismic  wave  fluctuations  at  the  surface  of  the  earth 
are 

a 

<«*>  =  2it*2J  dz  jdiKTW( 0,Kt.z)  sin’tKr 1  z /2k]  (7) 

ft 

<♦*>  -  2x*2|dz  \d1KT  W(0,Kt,z  )cos1(Jfj.2z/21t]  (8) 

t 

<u  ♦>  =  2rU2jdz  J d1Kr  W(0,Kr,z  ) 

srnltfr^/^lcoslArr^m]  (9) 


P 

Fig.  2.  Variance  ratio  y  versus  amplitude-phase  correlation  p  for  uni¬ 
form  extended  media.  The  solid  circles  are  from  Aki  [1973).  The  solid 
squares  are  from  Capon  [1974].  The  open  circle  comes  from  our 
analysis  of  NORSAR  data.  The  quantity  p  is  the  spectral  power  law 
index  (equation  (30)).  Position  along  one  of  the  curves  is  determined  by 
the  wave  parameter  D=4R  Ika  2  where  R  is  the  propagation  range,  k  is 
the  wave  manber,  and  a  is  the  correlation  length. 


An  example  of  information  to  be  gained  from  evaluation  of 
the  variances  is  given  in  Figure  2,  which  plots  the  two  quantities 


(10) 

(ID 


where  =  <u2>  and  =  <62>.  Both  y  and  p  are  ratios 
that  do  not  depend  on  the  overall  strength  <  p2  >  of  the  inho¬ 
mogeneities.  For  a  medium  that  extends  to  depth  R  with 
fluctuation  statistics  that  are  independent  of  depth,  both  y  and  p 
depend  on  one  parameter.  Given  a  scale  a  that  characterizes  the 
spectrum  W(  K)  ,  that  parameter  is  D  =  4R/(k  a1),  which  was 
defined  as  the  wave  parameter  by  Chernov  [I960].  For  each 
choice  of  medium  spectrum  there  is  a  unique  curve  in  y  -  p  space 
traced  out  through  different  values  of  D .  Some  examples  are 
given  in  Figure  2  along  with  some  experimental  results.  It  is 
seen  that  spectra  of  the  power  law  type  are  a  better  fit  to  the  vari¬ 
ances  than  a  Gaussian  spectrum  and  that  the  values  of  p  cluster 
around  0.3,  while  the  values  of  y  vary  substantially.  Note  that 
higher  seismic  frequencies  will  tend  to  have  lower  y  because  6m> 
increases  with  frequency  faster  than  u^.  Of  course,  a  medium 
whose  spectrum  is  dependent  in  some  way  on  depth  will  lead  to 
p,  y  values  that  do  not  lie  on  these  sample  curves,  but  the  general 
trend  of  smaller  p  for  power  law  spectra  versus  Gaussian  spectra 
will  apply. 

Further  information  beyond  the  variances  is  available  in 
coherence  functions.  For  example,  the  phase  coherence  between 
two  receivers  at  positions  x  ,r  and  xv  can  be  expressed  as 


<*(?.r  )*(?zr  )>  =  <<K?r)4>(0)>  (12) 


where  R  is  the  depth  of  the  deepest  inhomogeneities  being 
analyzed. 


where  we  have  defined  Xf  «  xjt  ■  xiT  and  we  have  used  hor¬ 
izontal  homogeneity  of  the  spectrum.  We  define  the  transverse 
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coherence  functions  (TCFs)  corresponding  to  each  of  the  above 
variances  <  u  ( xr )  u  (0)  >  ,  <  ♦  ( xT)  9  (0)  > .  and 
<  h  (  xT )  $  (0)  > ,  respectively.  Again,  since  horizontal  homo¬ 
geneity  of  the  spectrum  is  assumed,  these  functions  depend  only 
on  the  difference  between  the  two  positions.  The  expression  for 
the  TCFs  in  terms  of  the  medium  spectrum  are  (Flatte  and 
Moody,  manuscript  in  preparation,  1987) 

« 

<u(xT)u(0)>  =  2iU*|dzJd*Kr  W(0,ATr,z) 

sin*  cos( Kt  ■  xT )  (13) 

I 

<9(*r)9(0)>  =  Ink^dzfr  KTW{<S,KT,z) 

cos*[JCr*z/2*]cos[Kr  - xr]  (14) 

ft 

<«(xr)4>(0)>  =  2jtk*|  dz  jd1  KfW  (0,  KT,z) 

cos[JC,*z/2*]  sin[ATr*z/2Jt]  cos(Kr  •  ?r]  (15) 


Since  W  is  assumed  not  to  depend  on  the  direction  of  KT  but 
only  its  magnitude  (horizontal  isotropy),  the  KT  integral  can  be 
expressed  in  polar  coordinates,  and  the  angular  integral  may  be 
done.  The  results  are  then  a  function  only  of  die  magnitude  of 
the  separation  xT  and  not  its  direction: 

* 

<  u  (xT)  u  (0)  >  =  4it*k*|  dz 

tm\KT2iak]J<t,KTXT)  (16) 

a  - 

<<>(xr)4>(0)>  =  An2k}^dz^KTdKTW(Q,KT.z) 

cxx^Kj1  z  12k]  J a(KTXT)  (17) 

R  - 

<  u  (xr)$(0)>  =  4i^k^dz^KT  dKTW(0,KT.t) 

sm[Ar7-*r/2*]costArI-*r/2jtJ/0(XrJtr)  (18) 

We  find  it  convenient  to  work  with  normalized  TCFs  defined 
as  follows: 


KTdKTW(0,KT,z) 


<6(*r)9(0)>*  = 

<  ♦  C*T  )  ♦  (0)  > 

(19) 

<“  (xt)r  (0)>#  = 

<♦*> 

<**  (*r)“  (0)> 

(20) 

<“  (*r)6(0)>*  = 

<u  *> 

<  «  (*r)  u  (0)  > 

(21) 

These  normalized  functions,  like  die  quantities  y  and  p,  do  not 
depend  on  the  overall  strength  <  p*  >  of  the  medium  fluctua¬ 
tions.  We  note  that  the  values  of  these  coherences  it  xT  k  0  are 
1, 1,  and  p,  respectively. 

Note  that  our  expressions  for  the  variances  and  the  TCFs  are 
simpler  than  the  usual  Bom  approximation  volume-integration 
expressions  (S.M.  FI  site  and  T.  Moody,  manuscript  in  prepara¬ 
tion,  1987).  Here  we  have  a  one-dimensional  integral  along  z 
and  a  one-dimensional  integral  over  the  spectrum.  We  see  that 


i 


Fig.  3.  The  three  tnauvene  coherence  functions  (TCFs)  and  die  three 
angular  coherence  functions  (ACFi)  for  an  extended  medium  with  an 
exponential  correlation  function  with  icale  a  >10  km  and  a  aeiamic 
frequency  of  1  Hz.  The  different  curves  are  resulti  for  media  that  extend 
from  the  surface  to  different  depths  given  by  the  curve  labels  in  kQcme- 
ten:  that  i«,  varying  from  35  to  500  km. 


the  log  amplitude  fluctuation  u  is  weighted  toward  higher  wave 
number  KT  than  the  phase  fluctuation  due  to  the  "Fresnel  filter" 
factors  sin*tATr*z/2jk]  andcoc*[£rxzf2Jt].  We  also  see  that  this 
weighting  is  more  pronounced  for  smaller  r;  in  other  words, 
fluctuations  of  u  with  wavelengths  \g>  2x(tlk)H  are 
suppressed.  Note  that  there  is  more  suppression  near  die  sur¬ 
face;  it  is  difficult  for  wave  speed  fluctuations  to  cause  ampli¬ 
tude  fluctuations  without  a  sufficient  drift  space  before  the 
receiver. 

It  is  helpful  to  show  some  examples  of  the  normalized  TCFs 
for  an  interesting  medium  spectrum.  Figure  3  shows  the  three 
TCFs  for  a  medium  with  an  exponential  correlation  function 
with  scale  a  =  10  km.  That  is. 


W(K) 


<  p. *  >  a*  1 

it*  (1  +lf*o*)* 


(22) 


We  see  that  die  width  of  the  phase  coherence  function  is  larger 
than  the  width  of  the  medium  correlation  function,  while  the 
width  of  the  log  amplitude  coherence  function  is  smaller.  This 
is  due  to  the  Fresnel-filter  factor  which  favors  high  wave  number 
and  hence  small  scales  in  the  log  amplitude  coherence. 


3.  Angular  Correlation  Functions 

Consider  teleseismic  waves  arriving  from  two  different  direc¬ 
tions  indicated  by  unit  vectors  0i  and  9} ;  let  the  vectorial 
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Fig.  4.  Schematic  of  incident  P  waves  from  two  different  sources 
modem  an  a  seismic  array  (see  Figure  1). 


difference  in  angle  be  8  =  02  -  0,.  Finally,  the  magnitude  of 
the  angle  between  the  two  unit  vectors  is  denoted  by  9.  It  is  well 
known  that  the  patterns  of  u  and  p  on  the  Earth's  surface  will  be 
different  for  waves  bom  the  two  different  directions,  and  this 
difference  will  be  greater  for  larger  angular  differences:  that  is, 
larger  9.  (See  Figure  4.)  If  we  again  make  die  assumption  of 
horizontal  homogeneity  and  isotropy  of  the  medium  fluctuation 
spectrum,  then  statistical  correlation  between  the  log  amplitudes 
and  phases  of  die  waves  bom  two  different  directions  will  be  a 
function  only  of  the  magnitude  of  die  angular  difference 
between  the  two  waves;  that  is.  the  magnitude  of  9  and  not  its 
direction.  We  can  then  write  (S.M.  Flattd  and  T.  Moody, 
manuscript  in  preparation,  1987): 


<«  (8)«  (0)> 


a 


KTdKTW(0,KT,  z) 


]Jo(KTt  9) 


(23) 


<  4  (9)  0  (0)  >  =  4kV|  <fz|  Kt  dKr  W( 0, KT,z) 

ct»i[KT1zl'2k]Jo(KTz9)  (24) 

«• 

<  «  (9) 9  (0)  >  -  4****! dr| Kr  dKr  W( 0. Kr,z) 

sin(  ATr]r/2k  ]oos[  JCj-2z/2fc  ]J0(Krz  9)  (25) 


Thus  these  angular  correlation  functions  (ACFs)  are  very 
similar  to  dm  TCFs  except  that  the  transverse  separation  xT  has 
been  replaced  by  r  9;  in  other  words,  the  transverse  separation  is 
a  lines  function  of  depth,  given  by  the  separation  of  two  rays 


converging  to  the  same  point  on  the  surface  of  the  earth  bom  the 
two  different  directions  9j  and  02.  (See  Figure  4.) 

The  normalized  ACFs  have  at  least  two  important  properties: 
First,  they  provide  another  set  of  observables  bom  array  data 
that  give  a  different  weighting  to  inhomogeneities  at  different 
depths  than  the  TCFs  and  hence  provide  further  important  con¬ 
straints  on  medium  models.  Second,  the  variable  0  depends  on 
having  at  least  two  sources;  since  there  are  many  more  source 
(earthquake)  locations  than  there  are  stations  in  a  typical  seismic 
array  and  since  these  earthquakes  provide  angular  differences 
down  to  fractions  of  a  degree,  the  ACFs  provide  a  much  finer 
resolution  of  the  medium  inhomogeneities  than  the  TCFs.  For 
example,  the  NORSAR  subarrays  are  separated  by  distances  of 
10  km  or  more,  so  that  the  transverse  coherence  for  separations 
xT  smaller  than  about  10  km  cannot  be  observed.  However, 
observations  with  angular  differences  in  incoming  direction  of 
about  1°  are  not  uncommon,  providing  a  probe  of  scales  of  1  km 
even  at  a  depth  of  60  km. 

Predictions  for  normalized  ACFs  for  the  medium  spectrum  of 
(22)  (an  exponential  medium  correlation  function)  are  shown  in 
Figure  3.  Now  the  depth  to  which  the  medium  extends  becomes 
a  more  important  factor  in  determining  the  curves.  We  see  that 
the  ACFs  become  narrower  for  media  extending  to  larger  deprth 
if  other  parameters  are  kept  the  same. 

4.  Experimental  Results 

Previous  results  for  the  variances  have  been  given  by  A  fa 
[1973).  Capon  [1974J,  Berieussen  et  al.,  [1975a.  b],  Powell  and 
Mehzer  [1984],  and  others.  In  many  cases  they  presented  then- 
results  in  terms  of  a  uniform,  isotropric  random  layer  with  a 
Gaussian  sprectrum.  We  have  inverted  the  relevant  equations 
where  necessary  to  revive  their  variance  data  with  no  imposition 
of  a  medium  model.  In  cases  where  event-by-event  results  were 
given,  we  have  averaged  all  the  events  presented.  The  results 
are  given  in  Table  1.  We  see  general  agreement  between  the 
various  determinations  when  expressed  in  this  form,  except  for  a 
few  cases  where  a  very  small  number  of  events  were  used. 
Variance  results  are  graphed  in  Figure  2,  where  it  is  seen  that 
individual  events  have  a  reasonably  large  scatter,  averages  are 
consistent  with  a  power  law  typo  of  medium  rather  than  a  Gaus¬ 
sian  spectrum. 

Previous  results  for  TCFs  were  restricted  to  the  log  amplitude 
and  phase  without  the  cross-coherence  TCF.  Aid  [1973J  and 
Capon  [1974]  have  presented  approximate  TCFs  from  LASA 
data  showing  decorrelation  scales  of  the  order  of  10  km. 

Previous  results  for  ACFs  consist  only  of  qualitative  com¬ 
ments.  Berieussen  [1975]  remarks  that  events  are  well  corre¬ 
lated  in  amplitude  fluctuations  if  they  are  within  a  slowness 
difference  of  0.5  s/deg,  which  implies  that  their  incoming  direc- 

TABLE  I .  Experimental  Values  for  Variances 


Reference 

^naa* J 

unm 

P 

Aki  (1973) 

0.13 

0.32 

0.35 

Capon  [1974] 

0.10 

0.40 

- 

Berieussen  el  al.  [1975a,  6] 

0.006-0.11 

0.2-0.4 

- 

This  work 

0.135±0.004 

0.41040.006 

02640.03 
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E  -  W  Location  (km) 

Fig.  S.  Subarray  spatial  distribution  at  NORSAR.  Subamy  01 A  it  at 
latitude  60.82  N  and  longitude  10.83  E. 


lions  are  within  about  2°  and  gives  an  example  showing  that  two 
events  are  substantially  decorrelated  if  they  are  7°  apart  in 
azimuth,  which  implies  that  their  incoming  directions  are  5° 
apart. 

There  are  many  other  publications  that  include  attempts  to  fit 
a  deterministic  structure  to  the  data  from  one  or  many  events; 
we  will  not  try  to  discuss  their  data,  as  it  would  have  been 
presented  in  such  a  different  way  that  the  observations  relevant 
to  our  attempts  to  observe  small  scale  structure  would  not  have 
been  given.  We  discuss  (he  relationship  of  our  results  to  some 
of  the  deterministic  results  of  others  in  section  6. 

We  now  discuss  our  presentation  of  the  NORSAR  subarray 
data  given  in  two  NORSAR  reports  by  Berteussen  [1974]  and 
Berteussen  and  Husebye  [1974].  The  data  have  been  filtered  in 
two  important  ways.  First,  a  frequency  filter  has  restricted  the 
seismic  frequency  to  1-3  Hz,  so  the  nominal  center  frequency  is 
2  Hz.  Second,  the  elements  within  the  subarrays  have  been 
added  coherently;  each  subarray  has  about  six  elements  spread 
over  a  circle  of  radius  -33  km.  The  spatial  distribution  of 
subarrays  is  shown  in  Figure  5. 

The  travel  time  anomalies  given  by  Berteussen  [1974]  are 
given  for  104  beams;  each  beam  is  an  average  over  a  number  of 
events  whose  source  locations  are  close  such  that  their  incoming 
directions  at  the  NORSAR  array  are  within  a  few  tenths  of  a 
degree.  The  travel  time  anomalies  consist  of  arrival  time  resi¬ 
dues  with  respect  to  a  "best  plane  wave."  We  have  removed  the 
mean  arrival  time  for  each  beam  in  our  calculations.  Figure  6 
shows  the  arrival  angle  distribution  of  these  beams. 

The  subarray  intensity  variations  are  given  by  Berteussen  and 
Husebye  [1974]  for  183  beams.  We  have  normalized  these  data 
such  that  for  each  beam  the  mean  log  amplitude  is  zero. 

With  this  brief  description  of  the  gathering  of  the  data  (more 
details  are  given  in  the  original  reports)  we  now  describe  our 
method  of  estimating  variances,  TCFs  and  ACFs.  More  details, 
including  our  method  of  estimating  the  associated  statistical 
errors,  are  given  in  die  appendix. 

If  we  had  subarrays  that  were  far  apart  and  beams  that  were 
far  apart  in  angle,  then  the  information  from  each  subarray  and 


from  each  beam  would  be  statistically  independent  The  subar¬ 
rays  are  far  enough  apart  in  location  that  they  are  reasonably 
independent  for  the  purposes  of  error  calculations.  However, 
the  beams  are  often  close  together,  and  hence  the  assumption  of 
independent  beams  is  not  valid.  If  the  beams  are  not  indepen¬ 
dent  then  the  calculation  of  any  quantity  should  use  an 
appropriate  weight  for  each  beam.  Let  b  be  an  index  over  the  N 
beams;  then  we  define  w»  as  the  weight  of  beam  b .  The  appen¬ 
dix  describes  our  form  for  the  weight  In  order  to  calculate  rea¬ 
sonable  errors,  we  must  establish  the  effective  number  of 
independent  beams;  this  is  also  done  in  the  appendix.  We  note 
that  the  number  of  independent  beams  is  different  for  the  dif¬ 
ferent  TCFs;  we  find  42.1,  11.5,  and  32.4  for  the  effective 
number  of  independent  beams  in  determining  the  TCFs  of  log 
amplitude,  phase,  and  their  cross-TCF,  respectively.  The  ACFs 
have  similar  numbers,  although  the  particular  number  varies 
with  the  angle. 

To  describe  our  formulas  for  the  determination  of  experimen¬ 
tal  quantities,  consider  first  the  variances.  We  have 

N  kt 

<“J>  =  -xr  <26> 

<02>  -  77-  Z>*  (27) 

<“4»  =  ~  £  77  £  (28) 

where  b  is  an  index  over  the  N  beams,  Nw  =  and  s  is  an 
index  over  the  M  subarrays.  The  only  difficulty  in  the  above 
calculations  is  that  the  sets  of  beams  for  u  and  0  are  not  identi¬ 
cal.  In  order  to  calculate  the  cross  coherence,  we  have  searched 
through  the  two  sets  of  beams  and  found  corresponding  pairs 
that  are  within  1°  of  each  other  and  have  considered  them  the 
same  beam.  The  calculation  of  the  errors  is  somewhat  involved 
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Fig.  6.  Arrival  angle  distribution  of  the  various  beams  defined  at  NOR¬ 
SAR.  Each  beam  represents  data  from  a  duster  of  earthquakes  that  are 
close  enough  that  their  seismic  waves  arrive  at  NORSAR  with  angles 
that  are  within  a  few  tenths  of  a  degree.  The  circles  are  labelled  with  the 
angle  from  the  zenith  in  degrees.  The  axes  are  labelled  with  the  com¬ 
ponents  of  slowness. 


-19- 


6608 


Platts  and  W u-.  Smau,  Scale  Structure  Under  NORSAR 


0  20  40  80  80 

rT(hm) 


0  10  20  30 

6  (degrees) 


Pit-  7.  Observed  normalised  TCFi  and  ACFs  at  NORSAR.  Note  that 
the  regions  of  rapid  variation  in  the  TCFi  are  not  covered  well  by  the 
data,  became  the  different  subanays  have  a  minimum  separation  of 
about  15  km;  the  regions  of  rapid  variation  in  the  ACFi  are  well  covered 
by  the  data,  became  the  mmimtan  separation  in  angle  between  the  dif¬ 
ferent  beams  is  in  the  renge  of  a  few  tenths  of  a  degree. 


and  is  given  in  the  appendix.  The  results  of  the  calculations  of 
(26)-(28)  ate  given  in  Table  1. 

The  TCFs  are  calculated  with  formulas  analogous  to  (26)- 
(28).  In  order  to  plot  a  point  at  a  transverse  separation  xT,  we 
include  in  the  mm  all  pain  of  elements  with  their  separation  Xj 
within  a  given  bin.  The  remit*  of  this  procedure  are  shown  in 
Figure  7.  We  see  that  the  smallest  separation  available  is  about 
15  km,  and  the  coherence  is  already  below  05  at  that  separation. 

finally,  we  calculate  the  ACFs  by  the  tame  method,  with  the 
order  of  summation  reversed  and  with  selection  of  pain  of 
beams  whoa*  angular  difference  is  within  a  specified  bin.  The 
ACFs  are  also  shown  in  figure  7.  We  see  that  the  experimental 
points  in  die  ACFs  have  at  much  statistical  precision  and  more 
resolution  within  the  variable  renge  for  which  there  is  significant 
variation  in  das  coherence  than  the  TCFs.  The  higher  resolution 
of  the  ACFs  compared  with  die  TCFs  follows  from  the  closeness 
in  arrival  angle  of  the  beams  as  compared  to  the  spatial  separa¬ 
tion  of  the  tubarreys. 

Certain  dramatic  features  of  the  ACFs  provide  constraints  on 
earth  structure  models.  Tha  phase  coherence  function  has  a 
smooth  almost  linear  drop  to  a  small  coherence  at  an  angle  of 
about  13*.  The  log  amplitude  coherence  has  a  distinctively  dif¬ 
ferent  shape  dim  the  phase  coherence;  it  (hops  sharply  to  0.5  at 
an  angle  of  about  2*  but  retain*  significant  correlation  out  to 
angles  of  order  10*.  We  dull  see  in  the  next  section  the  impli¬ 
cations  of  these  features  for  Earth  models. 


It  is  to  be  emphasized  at  this  point  that  the  nature  of  the  NOR¬ 
SAR  array  restricts  the  scales  of  Earth  structure  that  can  be 
observed  with  NORSAR  data.  The  approximate  diameter  of  the 
array  is  1 10  km.  Because  means  and  "best  plane  waves"  have 
been  subtracted  from  the  data  across  the  array,  we  cannot 
observe  Earth  structure  with  scale*  larger  than  about  110  km. 
Because  the  data  consist  of  coherent  addition  of  stations  over 
subarrays  that  are  about  7  km  in  diameter,  we  cannot  observe 
structure  with  wavelengths  smaller  than  about  5  km.  Further¬ 
more,  our  typical  seismic  wavelength  at  2  Hz  of  4  km  also 
implies  that  we  cannot  observe  structure  wavelengths  less  than 
-4  km.  These  facts  have  to  be  taken  into  account  in  the  theoreti¬ 
cal  models  with  which  we  compare  the  experimental  data. 

It  is  true  that  the  separations  of  stations  within  subarrays  are 
in  the  1-  to  7 -km  range.  We  have  obtained  station  data  for  13 
nuclear  explosion  events  of  high  signal-to-noise  ratio  at  NOR¬ 
SAR.  We  have  found  that  the  station-by -station  time  anomalies 
are  less  than  0.1  s  and  are  difficult  to  measure  reliably  because 
of  the  differing  signal  temporal  structures.  The  signal  intensities 
are  more  reliable;  we  have  taken  a  5  s  time  window  at  the  start 
of  the  signal  and  have  calculated  the  log  amplitude  TCF  in  the 
same  manner  that  we  did  for  the  tubarray  data.  The  results  will 
be  shown  later  when  we  compare  models  to  the  data. 

It  is  of  interest  to  note  that  the  results  of  this  and  die  previous 
section,  that  is,  the  TCFs  and  the  ACFs,  can  be  used  to  under¬ 
stand  the  errors  and  possible  biases  of  sparse  arrays  in  measur¬ 
ing  die  magnitude  of  an  earthquake  or  nuclear  explosion.  The 
analysis  technique  of  McLaughlin  a  al.  [1988],  along  with  our 
experimentally  determined  coherence  functions,  could  be 
applied  to  this  problem. 

5.  Spectral  Models  Of  Earth  Structure 
Under  Nob sar 

We  picture  the  crust,  lithosphere,  and  asthenosphere  below 
the  Earth's  surface  as  having  random  variations  in  seismic  wave 
speed  that  are  described  by  a  spectrum  that  may  depend  on 
depth  and  may  be  anisotropic.  We  assume  that  the  spectrum 
does  not  depend  to  first  order  on  small  change*  in  geographical 
position  within  die  region  under  the  NORSAR  seismic  array. 
The  most  desirable  situation  would  be  one  in  which  we  had 
enough  experimental  information  to  invert  the  data  for  all  the 
characteristics  of  the  spectrum  describing  these  variations.  We 
do  not  delude  ourselves  into  thinking  that  our  present  data  set  is 
adequate  for  the  task.  However,  we  can  impose  some  con¬ 
straints  on  these  type  of  models  by  analysis  of  the  data,  and  we 
can  find  an  Earth  model  which  is  at  least  consistent  with  the  data 
and  with  what  we  know  erf  the  Earth's  structure  from  other 
information. 

We  need  to  modify  our  theoretical  expressions  for  the  fact 
that  the  NORSAR  array  has  a  finite  size  and  for  the  fact  that  the 
NORSAR  data  are  avenged  over  subarrays  and  is  at  finite 
nonzero  wavelength.  The  finite  size  is  modeled  by  removing 
from  the  spectrum  all  wavelengths  greater  than  1 10  km.  There 
is  a  considerable  uncertainty  in  exactly  where  this  upper 
wavelength  cutoff  should  be  because  it  depends  on  exactly  how 
die  "best  plane  waves”  were  removed  rod  how  averages  of  log 
amplitude  were  made;  the  main  comment  to  make  is  that  this 
cutoff  should  be  arranged  to  fit  foe  coherence  functions  at  large 
separations;  as  a  result  the  observed  large-separation  coherences 
are  not  a  strong  constraint  on  Earth  models.  This  cutoff  does  not 
introduce  artifacts  into  the  coherence  functions  at  small 
transverse  separations  or  at  retail  angular  separation*. 
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The  tubnny  averaging  might  be  modeled  by  assuming  that 
the  data  are  averaged  uniformly  over  the  area  of  a  circle  of 
radius  r„  =  3.5  1cm.  This  convolution  corresponds  in  the  spec¬ 
tral  domain  with  multiplying  our  spectra  by  the  following  func¬ 
tion  of  Kt  : 


C(ATr)e 


'VdKrr.)]2 

Krro 


(29) 


This  function  has  its  first  zero  at  KTrc  =  3.8,  or  a  medium- 
wavelength  scale  of  S.S  km.  In  our  theoretical  models  we  use  a 
small-wavelength  cutoff  of  5.5  km.  Again,  this  cutoff  does  not 
introduce  artifacts  into  the  coherence  functions,  but  it  does  affect 
their  exact  shape  at  very  small  separations  (of  order  5  km  or  1°); 
future  investigation  of  individual  station  waveforms  will  allow  a 
more  rigorous  treatment  of  these  cutoffs. 

The  first  model  we  may  try  is  one  with  an  exponential  correla¬ 
tion  function  with  a  scale  length  of  10  km,  following  those 
authors  who  have  looked  only  at  TCFs  [Ab',  1973;  Capon,  1974; 
Capon  and  Berteussen,  1974;  Berteussen  el  al.,  1975 a,  b; 
Powell  and  Meltter,  1984],  (See  Figure  8.)  We  see  that  the  TCF 
data  and  the  ACF  of  phase  are  not  badly  fit  with  medium  varia¬ 
tions  extending  down  to  about  200  km  but  the  ACF  of  log 
amplitude  and  the  cross  ACF  have  significant  disagreements 
between  model  and  data.  The  disagreements  have  to  do  with  the 
shapes  of  the  ACFs  in  the  regions  of  their  rapid  variations',  at 
angles  less  than  10°.  We  should  point  out  that  disagreements  at 
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Fig.  8.  Comparison  between  the  data  and  calculations  based  on  uniform 
random  extended  media  with  exponential  correlation  functions  having  a 
scale  length  of  10  km.  The  media  extend  from  the  surface  to  depths  of 
35  to  500  bn.  Note  that  the  best  medium  extent  for  the  ACF  of  phase  is 
about  200  km,  but  that  model  does  not  give  a  good  fit  to  the  shapes  of 
the  ACF  of  log  amplitude  and  the  cross  ACF. 
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Fig.  9.  Comparison  between  the  data  and  calculations  based  cm  uniform 
random  extended  media  with  power-law  spectra  having  index  p  as  3 
(self -similar).  The  media  extend  from  the  surface  to  depths  of  35  to  500 
km.  As  in  Fig.  8,  the  best  medium  extent  for  the  ACF  of  phase  is  about 
200  km.  This  self-similar  model  gives  a  better  fit  than  the  exponential 
correlation  model,  but  the  self -similar  model  still  does  not  fit  the  shape* 
of  the  ACF  of  log  amplitude  and  the  cron  ACF. 

large  angles  are  probably  less  important  because  the  shapes  of 
the  functions  there  are  affected  by  the  finite  aperture  of  the 
NORSAR  array.  It  is  of  interest  even  in  this  rather  unsatisfac¬ 
tory  model  that  structure  is  necessary  down  to  several  hundred 
kilometers  to  give  a  fit  with  reasonable  average  scale  lengths. 

The  data  in  Figure  2  indicate  (hat  a  power  law  spectrum  is 
also  appropriate  for  comparison  to  (he  data  (in  this  case,  y  and 
p).  We  will  therefore  restrict  further  models  to  a  particular  type 
of  spectrum;  a  power  law  expressed  in  the  form 

W(K)  =  AK->  (30) 

There  is  a  wide  variety  of  models  with  power  law  spectra 
because  the  power  law  index  p  may  change  with  depth,  as  may 
the  strength  A .  We  have  tried  many  examples  of  uniform  media 
down  to  a  cutoff  depth,  but  none  of  the  fits  were  judged  to  be 
acceptable.  For  example,  Frankel  and  Clayton  [1986]  suggested 
the  use  of  a  uniform  layer  with  a  power  law  medium  with  p  =  3, 
which  they  refer  to  as  a  self-similar  medium.  Figure  9  shows 
the  predictions  of  such  a  medium;  we  see  that  the  fit  is  better 
than  that  of  the  exponential  medium,  but  some  disagreements 
remain  in  the  ACF  of  log  amplitude  and  the  cross  ACF  at  angle! 
less  than  10°. 

We  have  found  one  combination  of  power  law  media  that  we 
feel  does  fit  the  data  reasonably  well.  We  have  not  yet  com¬ 
pleted  a  full  inversion  to  obtain  a  quantitative  idea  of  the  param- 
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Fig.  10.  Schematic  repreaentation  of  our  be«t  model  of  random  varia¬ 
tion*  under  the  NORSAR  array.  The  model  spectra  have  aharp  cutoffs  at 
low  and  high  wave  number  to  simulate  the  effects  of  the  finite  array 
aperture,  the  finite  size  of  the  subamyt,  and  the  finite  seismic 
wavelength. 


tier  uncertainties  in  this  model.  Our  best  model  consists  of  two 
overlapping  layers  (Figure  10).  The  first  layer  from  the  surface 
down  to  200-Ion  depth  has  a  power  law  index  p  =  0:  in  other 
words,  a  band-limited  white  spectrum.  The  second  layer  lie* 
between  depths  of  15  and  250  1cm  and  has  a  power  law  index 
p  =  4  and  a  strength  such  that  the  p  =  0  and  p  =  4  spectra  cross 
at  K  =  03 1  km~ 1 .  Note  that  p  =  4  corresponds  to  the  high  wave 
number  behavior  of  an  exponential  correlation  function. 

Comparison  of  the  data  and  our  best  model  is  shown  in  Figure 
11.  One  important  feature  of  the  data  that  constrains  the  model 
is  the  rapid  decrease  in  die  log  amplitude  ACF  compared  with 
the  phase  ACF.  Most  single-layer  models  have  a  similar  shape 
for  the  two  ACFs;  in  our  best  model,  the  log  amplitude  ACF  is 
controlled  by  the  shallow,  flat  spectrum  layer,  and  the  phase 
ACF  is  controlled  by  die  deep,  steep  spectrum  model. 

Figure  12  shows  the  predictions  of  several  variations  on  our 
beat  model:  the  shallow  layer  alone,  the  deep  layer  alone,  the 
result  of  extending  the  shallow  layer  to  250  km  instead  of  only 
200  km,  and  the  result  of  extending  the  deep  layer  up  to  the 
surface  instead  of  only  to  15-Ion  depth.  The  fits  are  clearly  infe¬ 
rior  to  the  beat  model,  although  die  differences  of  the  ltfter  two 
examples  from  our  best  fit  model  are  not  terribly  striking. 

The  rms  variations  in  phase  and  log  amplitude  are  directly 
proportional  to  the  rms  vviation  in  P  wave  velocity  in  our 
model.  However,  they  are  also  proportional  to  the  square  root  of 
the  medium  correlation  length  in  the  vertical,  which  we  have  not 
measured.  If  we  assume  that  the  medium  spectrum  is  isotropic, 
then  taking  the  observed  rms  value  of  u  as  0.41,  we  deduce  ■( 


rms  variation  in  P  wave  velocity  of  0.9%  in  the  upper  layer  and 
0.5%  in  the  lower  layer.  If  we  use  the  observed  rms  value  of  4. 
which  is  1.70,  we  find  values  of  2.2  and  1.3%.  We  regard  the 
difference  between  the  values  obtained  from  log  amplitude  vari¬ 
ance  and  from  arrival  time  variance  as  not  vary  significant 
because  the  determination  of  overall  variance,  as  opposed  to 
spectral  level,  is  notoriously  difficult  [ Flat! / et  at.,  1979].  The 
difference  does  imply  that  our  result  for  rms  P  wave  velocity 
fluctuation  must  be  stated  with  large  uncertainty;  that  is,  it  is 
between  1  and  4%,  if  an  isotropic  spectrum  is  assumed.  How¬ 
ever,  it  does  not  seem  likely  that  isotropy  is  appropriate.  If  we 
define  the  anisotropy  ratio  as  die  horizontal  correlation  length 
over  the  vertical  correlation  length,  then  our  result  for  the  velo¬ 
city  fluctuation  is  proportional  to  the  square  root  of  that  ratio. 
The  fact  that  the  upper  and  lower  layers  have  rms  variations 
within  a  factor  of  2  of  each  other  indicates  that  neither  are 
second-order  effects;  both  layers  are  significant  in  determining 
the  wave  fluctuations  observed  in  the  seismic  array. 

The  predictions  for  the  phase  and  cross  TCFs  fit  the  data  rea¬ 
sonably  well,  but  the  log  amplitude  TCF  prediction  appears, 
mainly  because  of  one  data  point,  to  decrease  more  rapidly  than 
the  data.  After  observing  this  discrepancy  we  used  13  nuclear 
explosion  events  at  NORSAR  to  determine  a  more  accurate 
experimental  log  amplitude  TCF.  (See  section  4.)  This  TCF 


Fig.  1 1.  Comparison  between  the  dau  aid  the  prediction!  of  oor  best 
model,  confuting  of  two  overlapping  Isyen.  Each  layer  hea  a  power  law 
apectrum  characterized  by  apectnl  index  p  .  The  character!  (tic*  of  the 
two  layer*  are  (0< 2  <  200  km;  p  =  0)  and  (15<Z  <250  km; 
p  »  4)  with  normalization  inch  that  in  the  overlapping  depth  region 
they  have  equal  apectnl  level*  at  wavenumber  0.31  km" The  ft  ia 
quite  reasonable,  particularly  in  the  region!  of  rapid  variation  Qeaa  than 
10*  in  the  ACFa). 
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12.  Companion  between  the  diu  and  the  predictioni  for  vsriitioni 
taken  from  our  beat  model.  Short  daihei:  the  deep  layer  only,  dotted: 
the  shallow  layer  only;  long  daihei:  our  belt  model  modified  by 
extending  the  ahallow  layer  to  a  depth  of  250  km;  aolid:  our  belt  model 
modified  by  extending  the  deep  layer  up  to  the  turface. 


involves  single  seismometer  stations  rather  than  subarray  aver¬ 
ages  and  hence  involves  separations  down  to  3  km.  We  were 
gratified  to  observe  the  agreement  between  our  best  model  pred¬ 
iction  and  this  data  set  obtained  after  we  had  developed  our  best 
model  mostly  from  die  ACFs.  (See  Figure  13.) 

6.  Geophysical  Interpretation  Of  The  Observed 
Medium  Variations 

The  Earth  is  characterized  by  vigorous  geological  activity 
driven  by  convection  in  the  mantle.  It  is  believed  that  most  of 
the  temperature  increase  with  depth  in  the  Earth  occurs  within 
thin  boundary  layers  at  the  top  and  bottom  of  the  mantle.  But 
this  is  a  statement  about  averages,  and  numerical  simulations  of 
mantle  convection  have  indicated  that  strong  variations  may 
occur  in  these  boundary  layers  [ Boss  and  Sacks,  198S;  Olson  et 
al  ,  1987].  These  variations  may  involve  temperature  differ¬ 
ences  of  up  to  1000  K,  with  spatial  scales  perhaps  comparable  to 
the  postulated  boundary  layer  thicknesses  of  order  100  km.  The 
variations  may  also  involve  compositional  differences  or,  in  the 
upper  mantle,  variations  due  to  partial  melting.  Although  there 
is  a  wide  variety  of  possibilities  for  the  formation  of  Earth  inho¬ 
mogeneities,  we  think  it  helpful  to  note  that  our  observations  can 
be  interpreted  as  small-scale  variations  driven  by  dynamic  man¬ 
tle  convection,  or  as  fossilized  compositional  differences, 
perhaps  induced  by  convective  processes  in  earlier  eras. 

The  change  from  the  shallow  layer  with  its  abundant  small- 


scale  structure  to  the  deep  layer  with  predominantly  larger-scale 
variations  could  be  indicative  of  some  change  in  rock  strength. 
One  way  to  obtain  small-scale  structure  is  from  a  statistical  dis¬ 
tribution  of  cracks,  but  it  seems  unlikely  that  cracks  would 
extend  to  1  OO-km  depths.  Another  way  would  be  to  have  a  ran¬ 
dom  distribution  of  relatively  homogeneous  layers  with  discon¬ 
tinuities  in  P  wave  speed  from  one  layer  to  the  next:  layers 
whose  thickness,  horizontal  extent,  and  wave  speed  discon¬ 
tinuity  across  interfaces  are  statistical  in  nature.  Our  results 
would  imply  that  these  random  layers  do  not  persist  below  about 
200-km  depth. 

The  rms  variation  in  seismic  wave  speed  that  is  implied  by 
our  spectral  model  is  a  few  percent  over  the  wavelength  band 
5.5-1 10  km,  which  is  obtained  from  the  variance  of  phase  or  log 
amplitude.  This  small  variation,  which  has  been  obtained  by 
previous  workers  as  well,  is  not  inconsistent  with  known  proper¬ 
ties  of  materials  combined  with  the  above  geophysical  parame¬ 
ters. 

We  may  note  some  relationships  of  our  results  to  previous 
studies  of  deterministic  structure  in  the  same  depth  regime.  All 
of  these  studies  deduced  structure  whose  lateral  variations  are  at 
scales  larger  than  100  km.  Our  data  discriminate  against  struc¬ 
ture  with  lateral  wavelengths  greater  than  110  km.  Therefore 
our  results  are  complementary  to  the  following  deterministic  stu¬ 
dies,  whose  results  should  in  some  cases  be  geologically  related 
to  ours. 

H addon  and  Husebye  [  1978]  constructed  a  model  with  a  sin¬ 
gle,  deep  dipping  layer  (thin  lens  model),  calculated  its  predic¬ 
tions  by  the  parabolic  wave  equation,  and  compared  them  with 
the  spatial  structure  of  arrival  time  and  amplitude  across  the 
NORSAR  array.  They  concluded  that  a  layer  at  a  depth  of 
150-200  km  explains  a  good  fraction  of  the  variances.  Their 


Fi*-  13.  The  TCF  of  log  amplitude  for  13  nuclear  explosion  events 
(eight  from  the  Soviet  Union,  five  from  the  United  States)  detected  at 
NORSAR  (open  circles).  The  data  were  calculated  from  individual- 
station  waveforms;  note  that  the  individual  stations  have  separations  that 
range  down  to  abcxx  3  km.  We  see  that  when  the  resolution  is  available, 
the  coherence  function  drops  to  less  than  0.5  in  a  few  kilometers.  For 
comparison,  the  solid  circles  are  the  same  points  shown  in  Figure  7, 
which  represent  subanay-averaged  data. 
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amplitude  data  had  been  imoothed  such  that  only  the  large-scale 
systematic  variation  was  kept.  Therefore  their  results  are  com¬ 
plementary  to  ours. 

Sacks  et  al  [1979]  deduce  a  discontinuity  at  230  Ion  from 
observations  of  long-period  precursors  to  direct  S  at  NORSAR. 
They  suggest  the  precursors  result  from  S  to  P  conversion  at  the 
discontinuity,  which  they  assume  to  be  the  lithosphere- 
astheno sphere  boundary.  Our  results  also  show  that  230  km  is 
an  interesting  depth.  We  would  rather  favor  the  interpretation 
that  230  km  is  die  bottom  of  a  thermal  boundary  layer  in  the 
mantle,  with  the  bottom  of  the  lithosphere  being  at  200  km. 

Given  and  Helmberger  [1980]  deduce  laterally  independent 
structure  of  P  wave  velocity  as  a  function  of  depth  for  northwest 
Eurasia,  using  short-period  and  long-period  body  waves  from 
nuclear  explosions.  They  see  a  low-velocity  zone  in  the  regime 
of  130-200  km,  which  is  likely  to  be  geologically  related  to  our 
observations  of  changes  at  200  and  230  km. 

Husebye  et  al.  [1986]  reported  a  seismic  tomographic  survey 
of  the  lithosphere-asthenosphere  beneath  southern  Scandinavia 
using  travel  time  residuals.  They  placed  the  boundary  between 
the  lithosphere  and  asthenosphere  between  100  and  200  km, 
with  lateral  variations  of  scale  greater  than  100  km. 

Thus  other  studies  of  the  larger -sc  ale  lateral  variations  under 
NORSAR  have  indicated  structure  at  depths  of  150-250  km, 
which  should  be  combined  with  our  observations  of  smaller- 
scale  lateral  variations  down  to  the  same  depth  range  in  order  to 
create  a  geologically  consistent  picture. 

Present  ability  to  model  the  detailed  dynamics  of  mantle  con¬ 
vection  is  so  crude  that  as  yet  no  meaningful  comparison 
between  the  strength  and  spectrum  of  variations  dial  we  observe 
and  the  expected  product  of  mantle  convection  is  possible. 
However,  one  can  hope  that  future  understanding  of  mantle  con¬ 
vection,  perhaps  involving  numerical  simulation,  will  be 
influenced  by  measurements  such  as  ours  of  die  strength  and 
spectrum  of  seismic  wave  speed  variations  in  the  thermal  boun¬ 
dary  layers  of  the  mantle.  More  generally,  improved  geological 
understanding  of  the  crust  and  lithosphere  may  come  form 
requiring  agreement  with  the  small-scale  variations  that  we 
observe  in  our  shallow  layer. 

Appendix 

TCFs  and  Their  Statistical  Errors 

We  begin  the  treatment  of  error  estimates  for  the  TCFs  by 
considering  the  expression  for  the  TCF  of  log  amplitude 
obtained  from  the  data  of  one  beam: 

c*(-*r)  *  XT  £  “*»“*<  (Al) 

where  p  is  an  index  over  the  M  pairs  of  subarrays  whose  separa¬ 
tion  is  within  the  desired  bin  of  xT.  (The  index  of  the  first 
member  of  the  pair  is  s ;  the  index  of  die  second  member  is  / .) 
The  approximate  variance  of  c»  is  obtained  as  follows:  First,  the 
expectation  value  of  c»  is  needed, 

<c»>  =  p,(xr)  (A2) 

where  p,  is  the  unnormalized  TCF  at  Xj.  Next,  the  expectation 
value  of  the  square  is  needed, 

<  ca2  >  “  - “w >  (A3) 

“pi 

Our  assumption  of  Gaussian  variables  allows  us  to  evaluate  the 
above  fourth  moment  in  terms  of  all  possible  permutation  of 


products  of  second  moments  [Flattd et  al.,  1979].  The  result  is 
then 

<e‘l>  =  P‘+  M<ul>1  (A4) 

where  we  have  ignored  many  of  the  off-diagonal  coherences  at 
small  compared  with  the  variance  of  h  .  Thus  the  variance  of  c* 
is 

(AS) 

It  is  important  to  note  that  because  we  have  many  beams,  we  can 
observe  at  each  Xf.  We  simply  plot  the  distribution  of  c* 
from  the  many  beams  and  calculate  the  variance  of  the  distribu¬ 
tion.  The  question  is,  what  is  the  uncertainty  in  ck  after  combin¬ 
ing  all  the  information  from  die  many  beams?  We  will  couch 
this  in  the  form  of  die  question  of  determining  N# ,  the  effective 
number  of  independent  beams. 

Since  we  have  data  from  N  beams,  we  should  combine  them 
in  an  appropriate  manner,  that  is,  with  weights: 

1  N 

c(xT)s  <  u  (xT)u  (0)>  =  —  £w»c»(xr)  (A6) 

Ww  *=1 

where  ft  it  an  index  over  the  N  beams,  is  the  weight 
assigned  to  each  beam,  and 

N.  =  (A7) 

»-i 

If  the  N  beams  were  independent,  then  the  optimal  combination 
would  require  wt  =  1,  and  it  is  simple  to  show  that  the  variance 
of  c  is  given  by 

al=ik<ul>1  (A8) 

We  define  the  number  of  effective  beams  as 

O*2 

(A9) 

o2 

That  is,  in  order  to  calculate  the  variance  of  our  measurement  of 
the  TCF  at  a  particular  xT,  we  observe  the  distribution  of  that 
TCF  from  our  many  beams,  determine  its  variance  from  that  dis¬ 
tribution,  and  then  divide  by  N,a.  For  the  case  of  independent 
beams,  from  (A3)  and  (A8)  we  see  that  is  equal  to  the  total 
number  of  beams  N . 

However,  many  of  the  beams  in  this  data  set  are  so  close  to 
each  other  that  the  information  in  each  is  not  independent. 
Therefore  a  more  accurate  measure  of  the  TCF  is  obtained  by 
weighting  each  beam  by  a  number  smaller  than  unity  if  there  are 
nearby  similar  beams.  We  have  chosen  weights  for  each  beam  in 
the  following  way: 

w»  = - ^ -  (A10) 

1  +  £p«(9w) 

irk 

where  p,  (0W )  is  the  ACF  between  beams  and  0*  is  die  angle 
between  beam  b  and  beam  d.  This  weighting  hat  the  appropri¬ 
ate  limits;  first,  it  is  equal  to  unity  if  all  beams  are  far  from  each 
other  and  therefore  independent;  second,  if  there  are  n  beams  on 
top  of  each  other,  each  beam  has  a  weight  of  1  In . 

Now  we  must  estimate  the  statistical  error  on  this  coherence 
function.  We  emphasize  that  our  treatment  of  errors  makes 
significant  approximations;  it  is  very  difficult  to  be  more  accu- 


-24- 


Flatte  and  Wu:  SmalvScalb  Structure  Under  NORSAR 


6613 


rate  in  the  evaluation  of  errors  without  going  through  an  exten¬ 
sive  program  of  numerical  simulation.  We  assume  that  all  the 
beams  are  in  a  narrow  cone,  so  that  the  expectation  value  of  the 
TCP  for  each  beam  is  the  same,  and  we  assume  that  the  statistics 
of  the  log  amplitude  or  arrival  time  are  close  enough  to  Gaussian 
to  estimate  the  error  reasonably.  Let  us  use  a  shorthand  notation 
in  which  xr  is  suppressed: 

1  H 

c  - —J^wtch  (All) 

*«1 

and  we  need  to  find  an  appropriate  expression  for  <c2>. 

1  N  N 

<c*>  =  r T-J  LEh'»h’s<cscs>  (A12) 

Nw  S-M-l 

But  each  ct  is  the  sum  of  products  of  two  log  amplitudes.  We 
may  indicate  this  by 

1  *  * 

<c^cd>  =  tty  ££<“  (*.  .0»  )u(y,  .04  )u  (xj  )u  (yy  ,6^  )>  (A13) 

where  Xjj  and  y,j  are  separated  in  space  by  Xj  and  0t  is  the 
angle  of  beam  b .  Thus  we  need  to  evaluate  a  fourth  moment  of 
log  amplitude.  The  Gaussian  assumption  provides  us  with  the 
means;  in  the  Gaussian  case  the  fourth  moment  is  the  sum  of  all 
possible  permutations  of  products  of  second  moments  { Flatt f  et 
al.,  1979].  Consider  first  the  permutation  in  which  the  pair  of  u 
from  beam  b  and  the  pair  of  u  from  beam  d  are  kept  together. 
Each  pair  has  an  expectation  value  of  p, ,  and  the  entire  contribu¬ 
tion  of  that  permutation  consists  of  p,2.  The  other  permutations 
involve  second  moments  with  one  u  from  beam  b  and  one  u 
from  beam  d.  We  consider  only  the  contributions  for  those 
pairs  of  u  that  are  al  the  same  station,  since  other  pairs  will  have 
a  substantially  smaller  coherence.  It  then  follows  that 

+  <AI4> 


%'  =  1  +  Ip.  (0^.(6!)  (A17) 

where  the  sum  is  over  all  the  other  pairs.  The  angles  9,  and  02 
are  chosen  so  that  the  product  of  the  p’a  is  the  maximum. 
Again,  this  form  has  the  correct  limiting  values. 

The  expression  for  the  ACF  can  be  given  then  as 

P.  (9)  =  tt-  T7  “«  (Al*) 

/v»  r-i  M  »-i 

where  q  and  r  are  the  indices  of  the  beams  making  up  pair  p 
and  s  is  the  index  of  the  station  within  the  array.  In  evaluating 
the  fourth  moment  as  a  product  of  second  moments,  we  again 
treat  the  first  permutation  completely  and  then  ignore  those  pairs 
which  belong  to  different  stations  in  subsequent  permutations. 
We  find 

Neff1  =  <A19) 

"»  p  i 

where  p  and  q  are  different  beam  pairs  with  angle  difference  0 
and  pM  is  the  normalized  ACF.  The  angle  0  is  the  smallest  angle 
obtained  by  selecting  one  beam  from  each  of  the  two  beam 
pairs.  The  angle  sy  is  the  angle  between  the  other  two  beams. 
Again,  this  is  an  approximation,  which  has  the  right  behavior  in 
the  limiting  case  of  many  beam  pairs,  each  of  which  is  made  up 
of  the  same  beams. 

It  is  of  interest  to  point  out  that  within  our  approximations,  the 
values  of  a  TCP  at  different  xT  and  the  values  of  an  ACF  at  dif¬ 
ferent  9  are  not  correlated. 
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_2  1  ££  P«2(»W) 

°  =  — I, —  <A15> 

'vw  ScU-1 


The  effective  number  of  beams  is  then  determined  from  (A5), 
(A9),  and  (A15)  to  be 


N«l  = 


N  N 


*=14=1 


pi  (0m) 
<«2>2 
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beams  instead  of  two. 
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